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2 Tom Irvine Background 


Tom Irvine earned B.S. and M.S. degrees in engineering 
science from Arizona State University, in 1985 and 1987, 
respectively. He wandered into the specialty of vibration 
serendipitously and found that vibration’s quirky traits were 
well-suited to his unconventional personality. 

He worked at Orbital Sciences Corporation in Chandler, 
Arizona, for 11 years as a full-time employee, and for 9 years 
as a consultant. He derived the shock and vibration levels and 
performed hands-on testing for the first Pegasus launch 
vehicle, which successfully flew in 1990. He has analyzed flight 
accelerometer data from Pegasus and numerous launch 
vehicles, as well as data from ground separation tests. 

He has posted hundreds of tutorial and software programs at his Vibrationdata website, covering shock, 
vibration, acoustics, structural dynamics, and signal processing. 

He later worked a scientist at Dynamic Concepts, Inc., in Huntsville, Alabama, where he performed 
vibroacoustic analysis on the NASA Space Launch System (SLS) Block IB Crew vehicle. He 
simultaneously served as an industry representative for the NASA Engineering & Safety Center (NESC) 
Loads and Dynamics Technical Disciplines Team. 

Tom is beginning a new vibrocoustic position at Blue Origin in Kent, Washington. 

Contact 

Tom can be contacted via Email: tom@irvinemail.org 

Readers are welcome to visit his blog: https://vibrationdata.wordpress.com/ 

Tom either performed or verified the calculations in this e-book using his Vibrationdata Matlab GUI 
package which is freely available at the above blog link. He hopes this book will help bridge the gap 
between shock and vibration analysis as taught in universities and as practiced in industry. 

Also note that this is expected to be the first of a series. Future works will include fatigue, finite 
element analysis, vibroacoustics, etc. 
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3 Variables 


A 

1. Amplitude 

2. Cross-sectional area 

A 

Acceleration 

B 

Bandwidth of an ideal rectangular filter 

C 

1. Damping coefficient matrix 

2. Crest factor 

C 

1. Damping coefficient 

2. Speed of sound in material 

c L 

Longitudinal wave speed 

Cp 

Bending wave phase speed 

c g 

Bending wave group speed 

D 

Plate stiffness factor 

D 

Diameter 

E 

Elastic modulus 

F 

Applied force in frequency domain 

F 

1. Frequency 

2. Applied force in time domain 

fn 

Natural Frequency (Hz) 

fr 

Ring frequency 

fs 

Strouhal frequency 

G,g 

Gravity acceleration, 9.81 m/sec 2 

H 

Transfer function 

H 

1. Impulse response function 

2. Plate thickness 

I 

1. Area moment of inertia 

2. Identity matrix 

k e 

Kinetic energy 

K 

1. Stiffness 

2. Wavenumber 

L 

Length 


M 

1. Mass matrix 

2. Bending moment 

m 

mass 

N 

Generalized displacement coordinate 
in frequency domain 

n 

1. Number of octaves 

2. Number of points 

3. Frequency index 

Pe 

Potential energy 

p 

Probability density function 

Q 

1. Amplification factor 

2. Eigenvector matrix 

q 

1. Eigenvector 

2. Generalized displacement 
coordinate in time domain 

s 

Strouhal number 

S 

Laplace variable 

T 

Period or duration 

t 

Time 

U 

Free stream air velocity 

U 

Longitudinal displacement 

V 

Shear force 

V 

Velocity 

w 

Wavelet or damped sine amplitude 

W 

Hanning window function 

X, x 

Displacement amplitude 

Y 

1. Amplitude 

2. Beam transverse displacement 

y 

Base displacement amplitude 

z 

Plate transverse displacement 
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z 

Relative displacement 

P 

Beam bending wavenumber 

5 

1. Logarithmic decrement 

2. Dirac Delta Function 

H 

Loss factor 

K 

Kurtosis 

X 

Wavelength 

M 

1. Poisson ratio 

2. Mean 

3 

Viscous damping ratio 

p 

Mass density (mass/volume) 

PL 

Mass density (mass/length) 

X 

1. Standard deviation 

2. Stress 

T 

Time delay in the convolution integral 

O 

Phase angle 

n 

Eigenvalue matrix 

n 

Angular frequency (rad/sec) 

cod 

Damped angular natural frequency 

(On 

Angular natural frequency 
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“If you want to find the 
secrets of the universe, think 
in terms of energy, 
frequency and vibration 

Nikola Tesla 
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4 Introduction 



Figure 4.1. IRVE 2 Launch & Free-Free Beam Analogy 


A rocket vehicle behaves as a free-free beam during flight. The vehicle’s body bending modes can be 
excited by wind gusts, aerodynamic buffeting, thrust offset, maneuvers, etc. The image shows the IRVE 2 
launch from Wallops Island. The vehicle is a Black Brant 9, which has a Terrier Mark 70 first stage. Flight 
accelerometer data showed that the fundamental bending frequency began at about 8 Hz and then swept 
up to 13 Hz as propellant mass was expelled. Autopilot guidance and control algorithms need to 
account for the body bending mode to maintain stability. 
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Figure 4.2. Tuning Fork, First Mode, A4 Note, Fundamental Frequency 440 Hz 

The mode shape from a finite element model is shown exaggerated. The two tines undergo an in-plane 
bending mode, 180 degrees out-of-phase with one another. The stem also participates in this mode, but 
its displacement is so relatively small that it is not apparent in the figure. 

4.1 Some History 

The Ancient Greek philosopher Pythagoras (570-495 BC) studied the vibration of stringed instruments 
and developed a theory of harmony. He first identified that the pitch of a musical note is in inverse 
proportion to the length of the string that produces it, and that intervals between harmonious sound 
frequencies form simple numerical ratios. 

The ancient Greeks believed that al objects were composed of some combination of the basic elements: 
earth, water, air, and fire. Aristotle (384-322 BC) attempted an explanation of earthquakes based on 
natural phenomena. He postulated that winds within the earth whipped up the occasional shaking of the 
earth's surface. He noted that earthquakes sometimes caused the water to burst forth in what would 
now be called a tsunami. 

Galileo Galilei (1564-1642) performed numerous experiments with oscillating pendulums during the 
Renaissance. He discovered that the period of swing of a pendulum was independent of its amplitude. 
He used his own pulse as a time measurement because there were no watches at that time. Christiaan 
Huygens (1629-1695) successfully built a pendulum clock based upon Galileo’s work. 
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Robert Hooke (1635-1703) developed the law of linear spring stiffness. This law has been generalized to 
the elasticity principle that strain in a body is proportional to the applied stress. 

Sir Isaac Newton (1643-1727) derived his laws of motion, which are the foundation of mechanical 
vibration analysis. He published his findings in a book known by its abbreviated title Principia. 

Jacob Bernoulli, Daniel Bernoulli, and Leonhard Euler derived the equation of motion for beam vibration 
circa 1750. 

Lord Rayleigh, John William Strutt, published Theory of Sound in two volumes during 1877-1878. 
Volume I covered harmonic vibrations, systems with one degree of freedom, vibrating systems in 
general, transverse vibrations of strings, longitudinal and torsional vibrations of bars, vibrations of 
membranes and plates, curved shells and plates, and electrical vibrations. Volume II covered aerial 
vibrations, vibrations in tubes, reflection and refraction of plane waves, general equations, theory of 
resonators, Laplace’s functions and acoustics, spherical sheets of air, vibration of solid bodies, and facts 
and theories of audition. 

The modern field of mechanical vibration analysis has been built upon the foundation of these authors’ 
works and has further developed as the result of failures and disasters, as well as the need to design, 
analyze and test structures and component to withstand dynamic environments. Vibration analysis is 
also important in other areas including semiconductor manufacturing, human exposure, energy 
harvesting, machine health monitoring, etc. 

4.2 Mechanical Failure 

4 . 2 .1 Failure & Damage Photo Gallery 



Figure 4.3. Turbine Blade Fatigue Failure 

A Qantas Boeing 747-400 aircraft was flying from Sydney to Singapore on May 9, 2011 when the crew 
noticed abnormalities from the aircrafts No. 4 engine during a climb. The indications included an 
increase in both the exhaust gas temperature and vibration levels. The plane continued to Singapore for 
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a safe landing and disembarkation of the passengers and crew. The Australian Transport Safety Bureau 
(ATSB) determined that the cause was a broken turbine blade. 



Figure 4.4. Washington Monument Crack. 


A one-inch wide, four-foot long crack formed in the Washington Monument, near the top of the 555- 
foot obelisk, due to the Mineral, Virginia earthquake on August 23, 2011. 

4.2.2 Historical Failures 

Railcar axles were failing under repeated “low level” cyclic stress, in the mid nineteenth century. These 
stresses puzzled engineers because the levels were much lower than the material yield stress. The 
Versailles rail accident occurred on May 8, 1842 after the leading locomotive broke an axle, and the 
carriages behind piled into it and caught fire causing many deaths. This prompted German scientist 
August Wohler to develop the S-N curves used in fatigue analysis. The term “fatigue” was chosen to 
describe the “tired metal” in the axles. 

The collapse of the Tacoma Narrows Bridge “Galloping Gertie” was captured on film on November 7, 
1940. This failure is often referred to as the classic “resonant vibration” failure. But it was more 
properly a “self-excited” or “flutter” response, as discussed in Section 22. 
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4.2.3 Aerospace and Other Industries 

Aerospace has myriad examples of potential vibration problems. Helicopters may undergo “ground 
resonance” prior to takeoff. Launch vehicles may have pogo oscillations in liquid engine propulsion 
systems. Solid rocket boosters may have thrust oscillations. Both high performance aircraft and launch 
vehicles must withstand random vibration due to turbulent boundary layers and shock waves as they 
accelerated through the transonic velocity and encounter the maximum dynamic pressure condition. 

Ships, automobiles, machine tools, buildings, nuclear reactors, and other mechanical systems and 
structures all have their own vibration concerns and failure modes. 

4.2.4 Failure Modes 

There are many types of potential failure modes including yielding, buckling, ultimate stress, fatigue, 
fretting, fastener loosening, relay chatter, and loss of sway space. Engineers must understand these 
hazards so that the components, systems, and structures may be designed and tested accordingly. 

Resonant excitation and is often a factor in these failures. Engineers thus have a responsibility to identify 
equipment natural frequencies through analysis and testing. The natural frequency is the frequency at 
which the system would oscillate if it were given an initial displacement and then allowed to vibrate 
freely. Resonance occurs when the excitation frequency is at or near the system’s natural frequency. 
Damping values, mode shapes, effective modal mass and other dynamic parameters are also needed for 
this analysis. 

4.3 Machine Health Monitoring 



Figure 4.5. Vibration Monitoring for Predictive Maintenance. 
(https://solutions.borderstates.com/ vibration-sensors-minimize-downtime/) 

Machines, pumps, and other equipment with reciprocating or rotating parts may experience vibration 
due to blade passing frequencies, gear mesh frequencies, magnetostriction motor hum, shaft 
misalignment, rotating imbalance, etc. In addition, fluid handling machines, like fans and pumps will 
experience broadband turbulence. 
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Some vibration is normal. But excessive vibration may cause accelerated wear and premature failure. 
High vibration levels could also indicate a bearing defect. Equipment manufacturers should provide 
acceptable limits in terms of amplitude and frequency. The amplitude specification may be in terms of 
displacement, velocity or acceleration. 

Vibration monitoring is thus needed to reduce maintenance costs, extend life, and improve safety. 
Machine vibration can be measured with accelerometers which are permanently mounted on the 
machine and monitored continuously with a wireless network. Or a technician may use a handheld 
device with an accelerometer that can be temporarily mounted to the machine via a magnetic base, as a 
periodic check. 

4.4 Vibration Energy Harvesting 

Piezoelectric transducers or electromagnetic induction devices can be used to convert ambient vibration 
energy to electrical current to charge batteries or to power wireless sensor networks. As an example, 
consider a cooling pump in a factory. The pump’s vibration could be used to power its own health 
monitoring sensor. 



Figure 4.6. Cantilever Piezoelectric Energy Harvesting Kit 
(https://www.mide.com/ products/ vlt-9001-piezoelectric-energy-harvesting-kit) 

Piezoelectric crystals are asymmetric. The asymmetry in the unit cell of the material sets up the 
mechanism whereby deforming the crystal leads to a small potential difference. The piezoelectric 
harvesters are typically cantilever beam structures. A mass may be added to the free end of the beam 
to tune the device to the source’s dominant vibration frequency. Note that steady vibration with a 
dominant frequency is best suited for harvesting. Outputs of 10 to 50 milliwatts are possible, depending 
on the vibration amplitude and frequency. 

Another potential benefit of an energy harvesting device is that it removes energy from a system, thus 
providing some damping which can extend the system’s life. 
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4.5 Human Vibration 

There are also biomedical concerns for the case where humans are exposed to vibration. Motion 
sickness and seasickness may result from vibration exposure in the 0.1 Hz to 0.5 Hz domain. This may 
occur because the brain has difficulty processing apparently conflicting data from the eyes, ear canals, 
and other somatosensory organs. 

Furthermore, operators of farm equipment, busses, and trains may suffer spinal damage due to long¬ 
term exposure. The human spinal column natural frequency is 10 to 12 Hz. Each organ and bodily part 
has its own natural frequency. Whole body vibration is addressed in ISO 2631 [I], Hand-arm vibration 
is another concern for operators of power tools. 

4.6 Vibration Excitation Sources 

4.6.1 Source Characterization 

Engineers must also determine the characteristics of the shock or vibration excitation so that they can 
properly analyze and test the affected components or structures. Design modifications can then be 
made to avoid dynamic coupling between the excitation frequencies and the structure’s natural 
frequencies. The excitation sources can be grouped into four types. 

4.6.2 Initial Displacement or Velocity 

A common example is a guitar string which is given an initial displacement and then suddenly released. 
Another is the Pegasus drop transient which is discussed in Section 10.6. An initial velocity example is a 
drop shock test machine. 

4.6.3 Applied Force 

Wind, acoustic pressure, and ocean currents impacting structures are all cases of applied force. So is 
footfall excitation on a pedestrian bridge. Shafts with rotating imbalance also impart an applied force 
within a machine. Another example is the force from a hammer striking an object. Impulse hammers 
are used to purposely excite a structure in modal testing to identify the structure’s natural frequencies, 
damping ratios and mode shapes. The force may also be applied by a small shaker attached to the 
structure via a stinger rod. 

4.6.4 Base Excitation 

Seismic excitation is the classic example of base excitation. Another is an automobile traveling over a 
speed bump or down a washboard road. Aircraft and launch vehicle avionics components are subjected 
to base excitation from their mounting structures during flight. These components are tested on large 
shaker tables to verify that they can withstand flight vibration prior to installation into the vehicle. 

4.6.5 Self-Excited Vibration 

Self-excited vibration is a special case where the alternating force that sustains the motion is created or 
controlled by the motion itself. This source type is also referred to as negative damping. Airfoil and 
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bridge flutter are two examples in this category. See the Tacoma Narrows Bridge failure in Section 22. 
The pogo oscillation in a launch vehicle with liquid engines is another case. 

4.6.6 Source Photo Gallery, Base Excitation 



Figure 4.7. Washboard Road, Base Excitation 
(https://www.trailerblocks.com/blogs/trailer-blocks-blog/tagged/trailer-tech) 

Ripples can form on gravel and dirt roads with dry, granular road materials. This pattern creates an 
uncomfortable ride for the occupants of traversing vehicles and hazardous driving conditions for vehicles 
that travel too fast to maintain traction and control. The resulting vibration may also damage 
suspension components. 

The earthquake caused the Cypress 
Viaduct to collapse, resulting in 42 
deaths. The Viaduct was a raised 
freeway which was part of the Nimitz 
freeway in Oakland, which is Interstate 
880. The Viaduct had two traffic decks. 
Resonant vibration caused 50 of the 
124 spans of the Viaduct to collapse. 
The reinforced concrete frames of 
those spans were mounted on weak 
soil. As a result, the natural frequency 
of those spans coincided with the 
excitation frequency of the earthquake 
ground motion. 

The Viaduct structure thus amplified the ground motion. The spans suffered increasing vertical motion. 
Cracks formed in the support frames. Finally, the upper roadway collapsed, slamming down on the 
lower road. The remaining spans which were mounted on firm soil withstood the earthquake. 



Figure 4.8. Loma Prieta Earthquake 1989, Base Excitation 
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Figure 4.9. Shaker Table Testing Lateral Base 
Excitation 

A small satellite is mounted to a slip table which 
is driven by a large electromagnetic shaker. The 
purpose of the test is to verify that the satellite 
can withstand the flight vibration that will be 
imparted by the launch vehicle. 



Figure 4.10. Shaker Table Testing, Vertical Base 
Excitation (courtesy Unholtz-Dickie) 

An equipment rack is mounted via wire rope 
isolators to an expander head which is mounted 
in turn on an electromagnetic shaker. 
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4.6.7 Source Photo Gallery, Applied Force 

A small shaker applies a force 
excitation to an automobile fender via 
a stinger rod. The applied force and 
resulting acceleration at the input 
point are measured by an impedance 
head transducer. Response 

accelerometers may be mounted at 
various locations on the vehicle. 



Figure 4.11. Automobile Modal Test, Applied Force 
(Image courtesy of the Modal Shop) 



Figure 4.12. Thumper Truck, Applied Force 


A thumper truck is a vehicle-mounted ground impact system which can be used to provide a seismic 
source to perform both reflection and refraction seismic surveys for oil, natural gas and mineral 
exploration. A heavy weight is raised by a hoist at the back of the truck and dropped about three 
meters, to impact the ground. The resulting ground waveforms are measured with geophones. Some 
thumpers use a technology called "Accelerated Weight Drop" (AWD), where high pressure gas is used 
to accelerate a heavy weight Hammer (5,000 kg) to hit a base plate coupled to the ground. 
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Figure 4.13. Space Shuttle Orbiter Ku-Band Antenna Dithering, Applied Force to Orbiter 


The Ku-band antenna is the disk below the lower Atlantis decal. Dither is a vibration employed in some 
mechanical systems to avoid stiction and to ensure smooth motion. Stiction is short for static friction. 
The antenna was dithered via a command signal at a frequency of 17 Hz to maintain its ability to 
smoothly search for Tracking and Data Relay Satellite System (TDRSS) satellites. There was a concern 
that the antenna’s resulting vibration would interfere with sensitive microgravity experiments such as 
crystal growth. 



Figure 4.14. Space Shuttle & Ares I Solid Rocket Boosters 
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Solid rocket boosters have elongated internal combustion cavities which act as organ pipes. Vortex 
shedding within the hot exhaust gases drives standing pressure waves inside these cavities. These 
pressure waves have a fundamental thrust oscillation frequency with integer harmonics. The cavities can 
be modeled as closed-closed pipes because the nozzle throat diameter is relatively small compared to 
the cavity diameter. The Space Shuttle boosters’ thrust oscillation frequency was 15 Hz. The Ares I 
vehicle was projected to have a 12 Hz oscillation, but this vehicle was cancelled prior to what would 
have been its first flight. Note that the speed of sound in the internal hot exhaust gas is about 3500 
ft/sec. 


4.6.8 Source Photo Gallery, Self-Excited Vibration 

Rocket vehicles with liquid engines, such as 
the Titan II, may experience combustion 
instability, which causes excessive vibration 
forces. This is a potential source of self- 
excited vibration, whereby the elastic vehicle 
structure and the propulsion system form a 
closed-loop feedback system. 

There are several types of combustion 
instability vibration effects. The most 
common effect is “Pogo,” which is similar to 
Pogo stick motion. In this case, a low 
frequency oscillation in the combustion 
chamber, or propellant feed system, excites 
the longitudinal vibration mode of the entire 
rocket vehicle, or some other structural 
mode. This may create a cyclical energy 
exchange between the vibration mode and 
the propulsion system oscillation. The 
problem may also be initiated when a wind 
gust or some other perturbation excites the 
vibration mode. This vibration in turn causes 
an oscillation in the propulsion system, which 
further excites the vibration. 

Placing accumulators in the fuel and oxidant 
lines to damp out the pressure fluctuations 
solves this Pogo problem. The accumulator 
contains a volume of gas that acts like a soft 
spring to reduce the propellant frequency to 
Figure 4.15. Titan II Gemini Vehicle, Pogo Oscillation well below that of critical structural 

frequencies. The accumulator volume must 
be carefully selected to meet this goal. 
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5 Single-Degree-of-Freedom Systems & Basic Concepts 

5.1 Free Vibration Model 



Figure 5.1. Single-degree-of-Freedom System, Spring-Mass Model 


The mass value is m. The viscous damping coefficient is c. The spring stiffness from Hooke’s law is k. 
The displacement is x. The velocity is x. The equation of motion is derived using Newton’s law. 


X F = mx 

(5.1) 

m x = —kx — cx 

(5.2) 

mx + cx + kx = 0 

(5.3) 

The resulting equation of motion is a second-order, ordinary differential equation, linear, homogenous 
with constant coefficients. The angular natural frequency in radians/sec is 

ii 

a 

3 

(5.4) 

The natural frequency in cycles/second or Hz is 


f n — C0 n j 271 

(5.5) 

The natural frequency is the frequency at which the system would oscillate if it 
displacement and then allowed to vibrate freely. The period T is the inverse of the 

were given an initial 
natural frequency. 

T — 1 / f n 

(5.6) 
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Furthermore, the damping coefficient divided by mass can be represented as 
c/m = 2^co n 

The corresponding amplification factor is 
Q = l/(2i;) 

The damped natural frequency is 

®d=®nV 1 -^ 2 


(5.7) 


(5.8) 


(5.9) 


The free vibration solution to equation (5.3) can found using a Laplace transform and applying the initial 
displacement and velocity per Reference [2]. The resulting displacement is 


sin(co d t)j 

The resulting time history is a damped sinusoid, similar to the Pegasus 
accelerometer data shown later in this document in Figure 10.7. 


x(t) = exp(-^© n t)l [x(0)]cos(co d t) 


+ 


x(0) + ^co n x(0) 




(5.10) 


drop transient flight 


5.2 Conservation of Energy Method 

The conservation of energy method is based on the time derivative of the total energy of an undamped 
system 


^(k e +p e ) = o 


(5.11) 


Apply this method to the SDOF system in Figure 5.1 to derive the equation of motion. 


d fl .9 1 


dt 


• 2 1 . 2 

—m x +—kx 


= 0 


mxx + kxx =0 


(5.12) 

(5.13) 


Dividing through by the velocity term yields the equation of motion for the undamped SDOF system. 


inx + kx = 0 


(5.14) 
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5.4 Rayleigh Method 

The Rayleigh method can be used to determine the fundamental frequency of an undamped system by 
setting the maximum kinetic energy Kg to the maximum potential energy Pg. 

[^E ] max - J max (5.1 5) 

The kinetic energy for the SDOF system in Figure 5.1 is 

K E=| m * 2 (5.16) 

The potential energy for the SDOF system is 



Assume a displacement of 

x(t) = Asin(co n t) 

The corresponding velocity is 

x(t) = co n Acos(co n t) 

The energy terms become 

K e = a> n 2 A 2 cos2 (®n0 

p E = \ kAl 

Substitute the two energy terms into equation (5.15). 

— m oo n 2 A 2 = —kA 2 
2 2 


(5.17) 

(5.18) 

(5.19) 

(5.20) 

(5.21) 

(5.22) 


Algebraic simplification yields the expected natural frequency 

®n =V k / m 


(5.23) 


5.5 Natural Frequency Example 

An avionics component is modeled as a solid mass per Figure 5.1. The goal is to mount the component 
via elastomeric isolator bushings so that the natural frequency is 30 Hz. The isolators will filter out high 
frequency vibration energy from the base excitation. Calculate the required isolator stiffness using 
equations (5.4) and (5.5). 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


22 



enDAQ 


Chapter: Single-Degree-of-Freedom Systems & Basic Concepts 


k = co n 2 m = (271 f n ) 2 m 


( 5 . 24 ) 


k = 


2n 


30 


cycles 


sec 


\2^| i iuf_2 • • ^ 


(51bm) 


llbf sec /in 
3861bm 


= 460 


lbf 


in 


( 5 . 25 ) 


The stiffness k in equation (5.25) is the total isolator stiffness. Now assume that the component will be 
mounted via four isolators in parallel. The individual isolator stiffness would then be 115 Ibf/in. The 
next step would be to determine whether isolators with that stiffness value are commercially available. 
Otherwise, some adjustment or compromise would be needed. 

5.6 Octave Rule 

A one octave increase in frequency means that the higher frequency is twice the lower frequency. 
Conversely, the lower frequency is one-half the higher one. The number of octaves n between any two 

frequencies f i and f 2 is 


n _ l°g(f 2 /f 1 ) 

log (2) 


( 5 . 26 ) 


Note that a piano keyboard has steps of one-twelfth octave counting both the black and white keys. 

The octave rule-of-thumb in mechanical vibration analysis states that there should be at least a one 
octave separation between two frequencies to mitigate dynamic coupling. For example, consider an 
SDOF system subjected to a harmonic force or base excitation. The system’s natural frequency should 
be tuned to less than one-half or greater than twice the excitation frequency. 

Isolation or low tuning would be the case where the system’s natural frequency was at most one-half 
the excitation frequency. High tuning occurs when the natural frequency is at least twice the excitation 
frequency. The system is said to be hard-mounted if the natural frequency significantly greater than the 
excitation frequency. 

As another example, consider two SDOF systems that are to be joined together. The respective natural 
frequencies should be separated by on octave prior to mounting one with other. The Pegasus launch 
vehicle in Section 10.6 has a natural frequency of about 10 Hz. The payload’s own natural frequency 
should be high-tuned to 20 Hz or more to reduce dynamic coupling effects. Note that tuning below 5 
Hz is a poor choice due to the possibility of high relative displacement of the payload within the vehicle 
fairing, as well as possible interference with the autopilot control algorithm stability. Further 
information about this frequency requirement can be found in the Pegasus Payload User Guide. 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


23 











enDAQ 


Chapter: Single-Degree-of-Freedom Systems & Basic Concepts 


5.7 Some Mathematics 

5.7.1 Complex Variables 

Complex variables with real and imaginary components are used extensively in signal analysis and 
structural dynamics. The imaginary component is expressed as a scale factor applied to j = y [—1 . 

Let C be a complex variable with real component A and imaginary component B. 

C = A + jB (5.27) 

The variable C is a vector in the complex, Euclidean plane. The magnitude is the norm. 

|c| = Va 2 +b 2 

The phase angle (j) is 

<|) = arctan(B / A) 


(5.28) 


(5.29) 


5.7.2 Euler’s Equation 

Euler’s equation is a complex exponential function used for Fourier transforms and for structural modal 
response. 


ex P(-j a ) = cos(a)± jsin(a) ( 5 .30) 

Complex exponentials can simplify trigonometry, because they are easier to manipulate than their 
sinusoidal components. 
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6 Two-Degree-of-Freedom Systems 

6.1.1 Grounded Two-DOF System, Equation of Motion 




k 2 (x j — x 2 ) -c 2 (xi-x 2 ) 




A 

m i 






— k I ^ 1 “C|X, 


Figure 6.1. Two-Degree-of-Freedom System, Spring-Mass Model 

Newton’s law can be applied to the system in Figure 6.1 to derive the equations of motion. The steps 
are omitted for brevity. The resulting pair of ordinary differential equations can be represented in 
matrix form and are coupled via the damping and stiffness matrices. 


~ m i 

0 

0 “ 

"*r 

+ 

C 1 +e 2 

~ c 2 

_ xr 

+ 

kj+k 2 

- k 2 _ 

" x i‘ 

_ 

1-1 

o o 

1 _ 1 

m 2 

_ x 2_ 


~ c 2 

c 2 _ 

>2_ 


~ k 2 

k 2 _ 

_ x 2_ 



( 6 . 1 ) 


A shorthand form is 


Mx + Cx + Kx = 0 


( 6 . 2 ) 
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The matrices and vectors are 


M = 

mj 0 

, c = 

~ c l+ c 2 

—c 2~ 

, K = 

ki +k 2 

“ k 2~ 

, X = 

~ x l~ 


_ 0 m 2 _ 


. ~ c 2 

c 2 . 


~ k 2 

k 2 . 


_ x 2_ 


(6.3) 


Consider the undamped, homogeneous form of equation (6.2). 


Mx + Kx = 0 


(6.4) 


6.1.2 Grounded Two-DOF System, Eigenvalues & Vectors 

Seek a harmonic solution of the form 


x = qexp(jcot) 

(6.5) 

x = jroqexp(jo)t) 

(6.6) 

x = -ro 2 qexp(jcot) 

(6.7) 


The q vector is the generalized coordinate vector. The undamped, homogeneous equation is 
transformed through substitution and algebraic manipulation into the generalized eigenvalue problem. 

jK-co 2 Mjq=0 (6.8) 

2 

The eigenvalues co can be found by setting the determinant equal to zero. 
det{K-w 2 Mj = 0 (6.9) 

There is an eigenvalue for each degree-of-freedom. Each angular natural frequency is then calculated 
from the square root of the respective eigenvalue. The corresponding eigenvectors represent 
orthogonal mode shapes. The eigenvector for each mode is found via 

{K-0) 2 Mjqi =0, i = 1,2 (6.10) 
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An eigenvector matrix Q can be formed. The eigenvectors are inserted in column format. 


Q = k <&] 


qi2 

c l21 ^22 


( 6 . 11 ) 


The coefficient q^j represents the modal displacement of mass i for mode j. 

Each eigenvector can be multiplied by an arbitrary scale factor. A mass-normalized eigenvector matrix 
Q can be calculated such that the following orthogonality relations are obtained. 


q t mq = i 

(6.12) 

q t kq = q 

(6.13) 


The superscript T represents matrix transpose. The identity matrix is I, a diagonal matrix of ones. The 
diagonal matrix of eigenvalues Q. 


6.1.3 Grounded Two-DOF System Example, Free Vibration 

Now define a modal coordinate rj(t) in terms of the normalized eigenvector matrix such that the 
displacement vector is 


x = Qr| = 


"qn 

j\2\ 


4l2 

^22 


Substitute equation (6.14) into equation (6.2). 


MQf) + CQf] + KQr| = 0 

Premultiply by the transpose of the normalized eigenvector matrix. 


Q T MQq + Q T CQf] + Q T KQrj = 0 
The orthogonality relationships yield 


(6.16) 


Ifj + Q T CQ f] + Qr| = 0 

Furthermore, the following simplifying assumption is made for the damping matrix. 


Q t cq= 


2^C0! 

0 


0 

2 ^ 2 w 2 


(6.14) 


(6.15) 


(6.17) 


(6.18) 
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The equation of motion can now be written as 


1 0“ 


4- 

"2^i coj 

0 

_ T|l" 

4- 

©! 2 

0 

_r ii" 


"0" 

0 1_ 

J\2_ 


0 

2 ^2 ®2. 

_ 1 l2_ 


0 

-1 

<N 

(N 

3 

_ r l2_ 


_0_ 


The two equations are now uncoupled in terms of the modal coordinates. The modal displacement for 
free vibration response to initial conditions is found via Laplace transforms as shown in Reference [3]. 
The modal displacement for dof i is 


T|i(t) = exp(-^j coj t )l rii (0)cos^ca d j t) + —cojrii (0) + f|i (0)]sin(co d i t) t (6.20) 

[ ®d,i J 

The damped angular natural frequency is 

®d,i = G W 1 -Si 2 ( 6 - 21 ) 

The physical displacement can then be found by substituting equation (6.20) into (6.14). The physical 
displacements can be expressed in terms of the eigenvector coefficients and modal coordinates as 

x l(t) = < lii T l 1 (t) + q 1 2 T l 2 (t) (6.22) 

x 2 (t) = q 2 i T| 1(0 + ^22*12(0 (6.23) 

The modal displacement initial conditions are required for the complete solution. These can be found 
from the following transformations. 

r[(0) = Q T Mx(0) ( 6 - 24 ) 

fj(0) = Q T Mi(0) ( 6 - 25 ) 


6.1.4 Grounded Two-DOF System, Participation Factors & Effective Modal Mass 

The participation factors and effective modal masses are indicators of how excitable the modes are given 
uniform base excitation to the given system. Some modes are more readily excited than others, and 
some cannot be excited at all for this excitation type. 
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The system’s generalized mass matrix m is given by 
m = Q T MQ 


(6.26) 


Again, the generalized mass will be the identity matrix if the eigenvectors are mass normalized. Let r 
be the influence vector which represents the displacements of the masses resulting from static 
application of a unit ground displacement. Define a coefficient vector L as 


L = Q T Mr 

The modal participation factor matrix Tj for mode i is 

r i Li / mu 


(6.27) 


(6.28) 


The effective modal mass m e ff ? x for mode i is 

m eff,i = Li 2 / A ii 


(6.29) 


6.1.5 Two-DOF System Example, Response to Initial Displacement 

The system in Figure 6.1 has the following parameters. 


Table 6.1. Grounded Two-DOF System, Parameters 


Variable 

Value 

mi 

2 Ibm 

m 2 

1 Ibm 

k i 

150 Ibf/in 

k 2 

100 Ibf/in 

X 1 

0.01 in 

x 2 

0.005 in 


Assume 5% damping for each mode. The initial velocity is zero for each mass. 

Apply the mass and stiffness values from Table 6.1 into equation (6.1) to form the following matrices. 



'2/386 

0 


' 250 

-100' 

M = 



, K = 




0 

1/386 


-100 

100 


(6.30) 
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The 1/386 factor is needed to convert Ibm to Ibf sec A 2/in. 

The eigenvalues and vectors are found by inserting these matrices into equations (6.9) and (6.10). The 
results are shown in the following table. 

Table 6.2. Grounded Two-DOF System Modal Parameters 


Mode 

Angular Natural 
Frequency (rad/sec) 

Natural 

Frequency (Hz) 

Participation 

Factor 

Modal Mass Ratio 

i 

125.25 

20.0 

0.085 

0.935 

2 

266.72 

42.5 

-0.022 

0.065 


The mass-normalized mode shapes in matrix format are 


8.93 -10.64 

15.05 12.63 


(6.31) 


The response time histories for the initial conditions are shown in Figure 6.2. 



Figure 6.2. Grounded, Two-DOF System, Response to Initial Displacement 
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6.1.6 Semi-definite, Two-DOF System 

6.1.6.1 Semi-definite, Two-DOF System Equation of Motion 


Xi 



-\A/W 



◄- k(Xl-X2) 





r 





k ( Xl - X2 ) * 

m 2 

Ungrounded, Two-DOF System 



X2 


Consider the system in the above figure. This system could represent a simple, two-node finite element 
model of a rod’s longitudinal vibration. A characteristic of this system is that the fundamental mode is a 
rigid-body mode at zero frequency. Both masses move in unison for the rigid-body mode. Newton’s 
law can be applied to the system to derive the equations of motion. The steps are omitted for brevity. 
The resulting pair of ordinary differential equations can be represented in matrix form. 


m. 

0 ~ 


+ 

' k 

-k" 

" x i" 


"O' 

0 

m 2 _ 

_ x 2_ 


-k 

k 

_ x 2_ 


_0_ 


(6.32) 


The eigenvalues and vectors are found using the method shown in Section 6.1.2. Each angular natural 
frequency is then calculated from the square root of the respective eigenvalue. 

G>1=0 (6.33) 


co 2 = 


' k(m 1 +m 2 ) 

m l m 2 


(6.34) 


The mass-normalized eigenvectors, representing the mode shapes, are 

1 


qi = 


V m 1 + m 2 


(6.35) 


q2 = 


_l_ 

Vm 1 m 2 (m 1 +m 2 ) 


m 2 
- m i 


(6.36) 
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Again, the two masses move in unison for the first mode. They move 180 degrees-out-phase for the 
second mode. The semi-definite system is revisited in Section 12.3.3. 


6.1.6.2 Semi-definite, Two-DOF System Example 

Consider the system in Figure 6.3 with the values in the following table. 


Table 6.3. Semi-definite, Two-DOF System, Parameters 


Variable 

Value 

mi 

2 Ibm 

m 2 

1 Ibm 

k 

2000 Ibf/in 


The mass and stiffness matrices are 



'2/386 

0 


' 2000 

-2000' 

M = 

0 

1/386 

, K = 

-2000 

2000 _ 


(6.37) 


The 1/386 factor is needed to convert Ibm to Ibf sec A 2/in. 

The eigenvalues and vectors are found by inserting these matrices into equations (6.9) and (6.10). The 
angular natural frequencies are 


(£>! = 0 

(6.38) 

0)2 =1076 rad/sec (171.3 Hz) 

(6.39) 

The mass-normalized eigenvectors in column format are 


, . H 1.343 -8.021' 


Q = [ qi 1 q? 1 = 

1 1 z J [11.343 16.042 

(6.40) 

The total mass sum is 3 Ibm (0.00777 Ibf sec A 2/in) which is also the effective modal 
body modes. The second mode has zero effective modal mass. 

mass of the rigid 
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7 Continuous Systems with Longitudinal Vibration 

7.1 Rod Longitudinal Natural Frequency 

7.1.1 One-Dimensional Longitudinal Vibration Equation of Motion 


E, A,p L 

u(x, t) 




“ L * 


Figure 7.1. Slender Rod 

The longitudinal displacement u(x,t) in a rod in undamped free vibration is governed by the second 
order, partial differential equation 


_d_ 

dx 


EA(x) 


8u 

dx 


= Pl( x ) 


d 2 u 

dt 2 


(7.1) 


The term EA(x) is the product of the elastic modulus and cross-sectional area. The equation for a 
uniform rod is 


d 2 u |f p^j d 2 u 
dx 2 2 


(7.2) 


Note that 


P = Pl/ a (7.3) 

Equation (7.2) is a common formula for one-dimensional wave propagation. Similar equations govern 
the propagation of sound in a pipe and the torsional vibration in a shaft. Note that a wave is the 
phenomenon in which physical energy propagates through space relative to a medium. Wave 
propagation is discussed throughout this document, as shown for seismic waves in Section 18.6.1 as for 
example. 
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The method of separation of variables can be applied to equation (7.2) as 

u(x,t) = U(x)exp(jcot) 

8 

—u(x,t) = j roU(x)exp( jcot) 



u(x,t) = -co 2 U(x)exp( jcot) 


(7.4) 

(7.5) 


(7.6) 


Substitute equations (7.4) and (7.6) into (7.2). The resulting spatial equation after some simplification is 


U(x) + co 2 

dx 2 




V w 


U(x) = 0 


(7.7) 


Longitudinal waves are non-dispersive. The group and phase wave speeds are equal. The wave speed is 
calculated from the elastic modulus and mass density. 


C L = V E /P 


The speed of sound is related to the frequency and wavelength by 


C T =f^ 


The spatial equation can be expressed as 


(7.8) 


(7.9) 


2 2 
d -U(X)+-^-r-U(x) =0 


dx" 


-L 


(7.10) 


The spatial solution is 


U(x) = ai sin 


r \ 

cox 


V C L j 


+ a2COS 


f \ 

cox 


v c l j 


(7.11) 


The a^ values in equation (7.4) are coefficients that depend on the boundary conditions. 
The boundary condition for a fixed left end is 

u(o) = o 


(7.12) 
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The boundary condition for a free left end is 


_d_ 

dx 



= 0 


(7.13) 


The same boundary conditions could be applied at the right end. The longitudinal natural frequencies are 
given in the following table. 


Table 7.1. Longitudinal Beam or Rod Natural Frequencies 


Configuration 

Natural Frequencies 

Mode Shapes 

Fixed-Free 

(2n -1 ) Ct 

n = l,2,3... 

f (2n-l)7rx^ 
U(x) = sin ---- 

2T 

v ZJ " J 

Free-Free 

f n =Y^. - = 1,2,3... 

U(x) = cos 

^nTix^ 

l L , 


Fixed-Fixed 

f n =!^S n = l,2,3... 

U(x) = sin 

^n7TX^ 

l L J 



The mode shapes are normalized to a value of one. Note that the Free-Free beam also has a rigid-body 
mode at zero frequency. Further solution details are given in Reference [4]. Also, longitudinal 
vibration is revisited in Section 19.2. 

The previous method of separation of variables is a modal solution approach. As an alternative, a wave 
solution can be applied to equation (7.2) as 


u(x,t) = Asin(kx-cot-(|)) 


(7.14) 


The wavenumber kis related to the angular frequency, longitudinal wave speed and wavelength as 


k = co/c = 2 tt/7, 


(7.15) 


7.1.2 Rod Longitudinal Vibration Energy Formulas 

Natural frequency and modes shapes can also be derived from a rod’s energy terms via the Rayleigh-Ritz 
or finite element method. The following formulas are given for reference. 
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The total strain or potential energy P E of a uniform rod is 

\2 

dx 

The total kinetric energy K E of a uniform rod is 

x2 


r L . 

f du l 

\ EA 

Jo 

IdxJ 


k t 


r L ^ A 


\ EA 

Jo 

IdxJ 


dx 


(7.16) 


(7.17) 


7.1.3 Beam Longitudinal Vibration Example 



The suspended beam can be idealized as a free-free beam for longitudinal vibration. The pendulum 
hammer is raised to some initial angular displacement and then released. The hammer strikes the 
beam’s end plate on the left end, delivering a force impulse. The longitudinal modes of the beam amplify 
the input energy and deliver a base input shock pulse to the test article mounted near the right end. 

Bai and Thatcher described this shock test method in Reference [5]. This method is still used today, 
although the pendulum hammer is often replaced by a pneumatic gun firing a projectile into the beam’s 
end plate. 

The goal was to choose the beam length so that its fundamental frequency would match the “knee 
frequency” of a shock response spectrum specification at 2000 Hz. (Further information on this type of 
specification is given in Section 18. The length for an aluminum beam to meet this goal is calculated as 


T 1 Ct 200,000 in / sec . _ 

L =-— = —- ! -- = 50 m « 4 ft 

2 f n 2 (2000 cycles / sec) 


(7.18) 


Note that the longitudinal wave speed in aluminum is about 200,000 in/sec. 
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7.2 Coil Spring Surge 



Figure 7.3. Automotive Valve Springs & Lifters 

(https://www.enginelabs.com/ engine-tech/ cam-valvetrain/when-do-i-know-my-valvetrain-is-in- 

di stress /) 

Spring surge is another form of longitudinal vibration, which arises because springs have both stiffness 
and inertia properties. This is in contrast to the usual simplifying assumption in vibration analysis that 
springs are massless. Springs used in high-speed machinery must have natural frequencies well in 
excessive of the frequency of motion that they control. Otherwise, the spring itself may resonate, 
resulting in loss of engine performance or even catastrophic failure. This can occur in automotive value 
springs when the engine camshaft RPM speed is increased well above the normal operating speed. Race 
car engines are at particular risk for this problem. Spring surge is also a potential concern when coil 
springs are used as isolation mounts for equipment. 

The following spring natural frequencies are taken from Reference [6]. Note the similarities with the 
longitudinal rod formulas in Table 7.1. 
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◄- L -► 



Figure 7.4. Fixed-Free Spring 


The fundamental frequency of the fixed-free spring is 


f n = 


2n-1 l 

1 

4 J 

V 


PlL 


, n = 1,2,3... 


(7.19) 


Its corresponding normalized mode shape is 
^(2n — l)7tx^ 


U(x) = sin 


2L 


(7.20) 


J 


The wave speed in the spring is 

Cl = V kL /pL 

The wave speed equation applies for other spring boundary conditions as well. 


(7.21) 


< - L-► 



I- * 

Figure 7.5. Fixed-Fixed Spring 
The fundamental frequency of the fixed-fixed spring is 




n = 1,2,3... 


(7.22) 
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The corresponding normalized mode shape is 
U(x) — sin 


^nTix^ 


J- J 


(7.23) 


Note that the Fixed-Fixed and Free-Free cases have the same natural frequency equation. 



Figure 7.6. Spring with End Mass 

The simplified fundamental frequency formula for a spring with and end mass M is 


fn = 


1 


27t ]l ^4 ■+■ 0.33p L L 
A more rigorous formula is 


The corresponding mode shape is 

sin(7. n x/L), n = 1,2,3,.. 


M 

PlC 


A, n , n = 1,2,3. 


(7.24) 


(7.25) 


(7.26) 
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7.3 Ring Frequency 



Undeformed Ring 


Radial Displacement, 
all points in-phase 


Figure 7.7 Ring Mode from Finite Element Model 


Consider a thin ring with a rectangular cross section and with completely free boundary conditions. The 
ring frequency corresponds to the mode in which all points move radially outward together and then 
radially inward together. This is the first extension mode. It is analogous to a longitudinal mode in a 
rod. The ring frequency is the frequency at which the longitudinal wavelength in the skin material is 
equal to the vehicle circumference. 


fr 



(7.27) 


The ring frequency is an important concept for launch vehicles. The front end of a launch vehicle is 
typically composed of cylindrical shell segments which house avionics components and the payload. 
These shells tend to have notable vibration responses at their respective ring frequencies due to 
external acoustic environments and stage separation shock events. Furthermore, there is a rule-of- 
thumb in statistical energy analysis (SEA) that a cylinder tends to behave as a flat plate above its ring 
frequencies. This is a simplifying assumption. Also note that the ring frequency is an idealized concept 
for a cylindrical shell. In practice, cylindrical shells tend to have a high modal density near the ring 
frequency. 

Consider an aluminum launch vehicle cylindrical shell with a diameter of 72 inches, 
longitudinal wave speed in aluminum is approximately 200,000 in/sec. The ring frequency is 

f r = (200,000 in/sec)/(7t60in) = 1060Hz 


Again, the 

(7.28) 
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8 Bernoulli-Euler Beams 


8.1.1 Beam Equations 


I 


El, pL 


y(x, t) 


Figure 8.1. Cantilever or Fixed-Fixed Beam 


El, pl 

hr 


I. 



Figure 8.2. Simply-Supported or Pinned-Pinned Beam 
The governing equation for beam bending free vibration is a fourth order, partial differential equation. 


d 2 y 


dx / 


EI(x) 


dx 2 


= Pl 


d^y 

dt 2 


( 8 . 1 ) 


The term EI(x) is the stiffness which is the product of the elastic modulus and area moment of inertia. 
The equation for a uniform beam is 


-El 




( 8 . 2 ) 


The method of separation of variables can be applied as 

y(x,t) = Y(x)exp(jcot) 


(8.3) 
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dt 


y(x,t) = jwY(x)exp(jcot) 


(8.4) 


dv 


-y(x,t) = -oo 2 Y(x)exp(jcot) 


(8.5) 


Substitution into equation (8.2) leads to spatial ordinary differential equations 


d -Y(x)-® 2 ]-^Y(x) = 0 


cbr 


El 


( 8 . 6 ) 


Define the wavenumber p as 




(8.7) 


The differential equation can be expressed as 


i4 

-Y(x)-p 4 Y(x) = 0 


dx 


( 8 . 8 ) 


The spatial solution has the form 


Y(x) = ai sinh(px) + a 2 cosh(Px) + a 3 sin(Px) + a 4 cos(Px) 


(8.9) 


Substitute the spatial solution into the differential equation. 


P 4 |aj sinh(Px) + a 2 cosh(Px) + a 3 sin(Px) + a 4 cos(Px)} 


co 2 j — Uaj sinh(Px) + a 2 cosh(Px) + a 3 sin(Px) + a 4 cos(px)} = 0 


El 


( 8 . 10 ) 


Equation (8.10) is satisfied by the wavenumber relationship in equation (8.7), which gives credibility to 
the assume spatial solution. The aj values in equations (8.9) and (8.10) are coefficients that depend on 

the boundary conditions discussed in 8.1.2. The spatial solution gives eigen function modes shapes. It 
also gives eigenvalue roots from which the natural frequencies are calculated. 
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The eigenvalues are represented as P n L . The angular natural frequencies C0 n are 

®n =Pn 2 V EI /pL 


( 8 . 11 ) 


8.1.2 Common Boundary Conditions 
8.1.2.1 Cantilever Beam 

Recall the cantilever or fixed-free beam in Figure 8.1. The displacement and slope at the fixed end are 
both zero. 


Y(0) = 0 


—Y(0) = 0 
dx 


( 8 . 12 ) 


The moment and shear force at the free end are both zero. 


M: 


EI(x)—rY(x) 
dx 2 


x=L 


v=— 

dx 


EI(x)—rY(x) 
dx 2 


(8.13) 


x=L 


Application of these boundary conditions to the spatial displacement in equation (8.9) and its derivatives 
yields the following equation for finding the coefficients. 


sin(pL) + sinh(pL) 

cos (PL) + cosh (pL) 

_a i" 


'o' 

cos (PL) + cosh (pL) 

-sin(pL) + sinh(pL) 

_ a 2_ 


0 


(8.14) 


By inspection, equation (8.14) can only be satisfied if = 0 and a 2 =0. Set the determinant of the 

matrix to zero in order to obtain a nontrivial solution. This yields the following transcendental equation 
for finding the roots by the Newton-Raphson method. There are multiple roots which satisfy equation. 
A subscript is thus added. 

cos(P n L)cosh(P n L) = -l (8.15) 

The fundamental root for the cantilever beam is |3 n L = 1.87510. The natural frequencies can be found 
by substituting each root in to equation (8.11). 
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8.1.2.2 Pinned-Pinned Beam 

Recall the pinned-pinned beam in Figure 8.2. The boundary conditions at the left end of the pinned- 
pinned beam are 


= 0 (8.16) 

x=0 

The boundary conditions at the right end are similar. 

Application of these boundary conditions to the spatial displacement in equation (8.9) and its derivatives 
yields the following equation for finding the coefficients. 


Y(0) = 0 , M = 


EI(x)—Y(x) 
dx 2 


sinh(pL) 

sin(pL) 

~ a l~ 


'O' 

sinh(pL) 

-sin(pL) 

_ a 2_ 


0 


(8.17) 


By inspection, equation (8.14) can only be satisfied if aj =0 and a 2 =0. Set the determinant of the 

matrix to zero in order to obtain a nontrivial solution. This yields the following transcendental 
equation. There are multiple roots which satisfy equation. A subscript is thus added. 

sin(P n L)sinh(p n L) = 0 (8.18) 

The equation is satisfied if 

P n L = rm, n =1,2,3,.... ( 8 | 9 ) 

The natural frequencies can be found by substituting each root in to equation (8.11). 

8.1.3 Beam Bending Fundamental Frequencies 

The solutions to the beam equation for various boundary conditions are derived in Reference [7]. A 
summary of beam formulas is given in Table 8.1. Note that the free-free beam also has a rigid-body 
mode with zero frequency. 
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Table 8.1 Beam Bending Fundamental Frequencies 


Configuration 

Fundamental Frequency 

Cantilever (Fixed-Free) 

i 

2n 

'3.5156" 
. L 2 

V ei /Pl 

Cantilever with end mass 



1 

1 3EI 

2n' 

V (0.2235 m beam + m en( j) L 

Pinned-Pinned 

1 

27t 

r -\ 2 

T « 

Free-Free & Fixed-Fixed 

1 

2n 

'22.373“ 

L 2 

V ei /Pl 

Fixed-Pinned 

1 

27t 

1 1 

00 

v5 J 

Isi 

l-H 

"O 

r 


The pinned-pinned beam has integer harmonics as follows. 


2n 



n = 1,2,3... 


( 8 . 20 ) 


The beam bending frequencies for other configurations have non-integer harmonics. The free-free beam 
higher modal frequencies are shown in Table 8.2. The free-free formulas also apply to the fixed-fixed 
beam. The fixed-free beam higher modal frequencies are shown in the following table. 


Table 8.2 Free-Free & Fixed-Fixed Beam Bending Frequencies, Elastic Modes 


Mode 

Natural Frequency (Hz) 

ft 

1 r 22.3731 /If 

2tt | l 2 J V P 

f 2 

= 2.757f i 

f 3 

= 5.404 f ! 

u 

= 8.402 f i 
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Table 8.3 Fixed-Free Beam Bending Frequencies, Elastic Modes 


Mode 

Natural Frequency (Hz) 

fl 

1 3.5156 r—- 

~ 2rc| L 2 

f 2 

= 6.268 f x 

f 3 

= 17.456 f i 

u 

= 34.386 fj 


8.1.4 Beam Bending Participation Factors & Effective Modal Mass 

Again, the mode shapes Y n (x) and their corresponding mode shapes are found by applying the boundary 
conditions to the displacement shape in equation (8.9). The participation factors are then 


r n =J 0 L pL Y n( x )dx 


The effective modal mass is 


m eff, n 


r L 

J 0 Pl Y n( x )dx 
Jq" Pl [ Y n( x )] 2 dx 


( 8 . 21 ) 


( 8 . 22 ) 
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8.1.5 Bending Wave Speed & Wavelength 



Distance 

Figure 8.3. Traveling Bending Wave Packet Dispersion, Snapshots at Three Different Times 

Both the red and blue circles oscillate in the vertical axis in the above figure. The red circle travels at 
the phase velocity. The blue circle propagates twice as fast at the group velocity and is always at or near 
the positive or negative peak in the wave packet. The phase speed is the more important of the two 
speeds for vibroacoustics analyses where the bending waves are excited by the external sound field, or 
vice versa. 

Bending waves are dispersive as a result of the governing fourth-order partial differential equation. The 
wave speed varies with the angular natural frequency. The phase speed for a given frequency is 

c p =V® [EI/p l ] 1/4 (8.23) 

The group speed is twice the phase speed. 
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c g“ 2c P (8.24) 

Demonstrating bending wave propagation in terms of these two speeds is best done with an animation, 
but the still images in Figure 8.3 may be useful. Furthermore, the phase speed is related to the frequency 
and bending wavelength by 

c p =fk (8.25) 


8.1.6 Beam Bending Energy Formulas 

Natural frequency and modes shapes can also be derived from a beam’s energy terms via the Rayleigh- 
Ritz or finite element method. The following formulas are given for reference. The total strain or 
potential energy P E of a uniform beam is 


c L 

f El 

fd 2y l 

O 

—» 

ldx 2 J 


dx 


The total kinetric energy K E of a uniform beam is 
Ke =^®„ 2 J 0 L p[y] 2 dx 


(8.26) 


(8.27) 


8.1.7 Beam Example, Wind Chimes 



The wind chimes in the figure are made from anodized aluminum. 
The tallest one is 34.7 inch long, 1.25 inch diameter, and 0.094 inch 
wall thickness. The chimes are like free-free beams undergoing 
bending vibration when struck. The vibration then generates 
sound waves at discrete frequencies matching the bending 
frequencies. The fundamental frequency of the tallest chime is 
calculated via the following steps. 


Figure 8.4. Author’s Wind Chimes 
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The area moment of intertia is 


i=-^ 

64 L 


d 0 4 -di 4 


64 n 


(1.25 in) 4 -(1.16 in) 4 


= 0.0321in 4 


(8.28) 


The mass per length is 


PL = 


o.i ^ 

V in j 


1 lbf sec z / in 
386 lbm 




v^y 


(l.25in) 2 -(l,16in) 2 


= 4.59e-05 


lbf sec 


in 


(8.29) 


The natural frequency is 


fn = 


2n 


22373 
(34.7in) 2 


lfl.0e + 07 ^ 

V in 1 J 


(o.0321in 4 ) j 4. 


59e-05 


lbf sec' 


= 247 Hz 


in 


(8.30) 


The tallest wind chime was struck separately from the others. The resulting sound was recorded and 
analyzed. The sound tone represented bending modes. The theoretical value from equation (8.30) 
agrees closely with the measured frequency 244 Hz. A comparison of the first four bending modes and 
their respective theoretical values from Table 8.2 is given in Table 8.4. There are two orthogonal planes 
where transverse vibration occurs. Each mode actually represents a pair of modes. The pair should 
ideally have the same frequency but may be slightly different due to manufacturing and material 
imperfections, as well as two holes in each for the connecting cord. 


Table 8.4. Tallest Chime Elastic Mode Natural Frequencies 


Mode 

Calculated (Hz) 

Measured (Hz) 

Musical Note 

1 

238 

244 

B 

2 

657 

664 

E 

3 

1288 

1273 

D# 

4 

2127 

2050 

C 
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Figure 8.5. Wind Chime Elastic Mode Shapes 

The Y-axis represents transverse displacement which has been normalized to have a maximum absolute 
value of one. The points where each curve crosses the zero baseline are nodal points. 
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Figure 8.6. Tallest Wind Chime Time History 

The time history shows some beat frequency effect because each mode is really a pair of orthogonal 
modes with closely-spaced frequencies. 
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Figure 8.10. Tallest Wind Chime Fourier Transform 

The highest response is at the fundamental bending mode at 244 Hz. The fourth bending modes at 2050 
Hz is barely excited. Note that Fourier transforms are covered in Section 15. 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


52 










enDAQ 


Chapter: Thin Plate Vibration 


9 Thin Plate Vibration 

9.1.1 Rectangular Plate, Bending Vibration 

9.1.1.1 Rectangular Plate Equations 



Figure 9.1. Rectangular Plate 

The fourth-order governing equation of motion for the out-of-plane displacement is 


D 


^4 
8 z 


+ 2 - 


,4 a 4 ^ 

8 z 8 z 


dx 2 dy 2 3y' 


. d 2 z 

+ ph—r- = 0 

at 2 


(9.1) 


The plate stiffness factor is 


D = 


Eh J 


12(l — pi 2 ) 


(9.2) 


The kinetic energy in the plate is 


Kr — 


phfr r b ra 2 


Uo z2<Jxdy 


2 JO JO 


(9.3) 
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The potential or strain energy in the plate is 


p _ D f b f a 

P e-yJ 0 Jo 


/ 9 \ 

d 2 Z 


dx z 


+ 


C o \ 

a 2 z 


a y 2 


+ 2p 


(TZ 

dx 2 


(TZ 

Sy 2 , 


+ 


2(1-M) 


a 2 z x 

8x8y 


dxdy 


(9.4) 


Assume a harmonic displacement response for the partial differential equation (9.1). 


z(x,y,t) = Z(x,y)exp(jcot) 


(9.5) 


Substitution of the assumed response leads to the spatial equation 


D 


r d 4 Z d 4 Z 9 4 Z^ 

-r + 2---- +-- 


dx* dx 2 dy 2 dy 


- phoc> Z = 0 


(9.6) 


J 


The modal solution equation (9.6) for most boundary conditions requires an approximation method 
such as the Rayleigh-Ritz or finite element method. Further information on plate vibration solutions is 
given in Reference [8]. 

9.1.1.2 Rectangular Plate, Simply-Supported along All Edges 

The case of plate with simply-supported boundary conditions on all four edges can be solved for exactly. 
The boundary conditions for this plate are 


Z(x,y) 

= 0, 

M x (x,y) = 0 

for x = 0, a 

(9.7) 

Z(x,y) 

= 0, 

M y (x, y) = 0 

for y = 0, b 

(9.8) 


The displacement and bending moment along each edge is equal to zero. 

Assume the following displacement function which satisfies the boundary conditions, where mand nare 
indices. The mass-normalized modes shapes for the plate are 


7 - 

^mn 


yj pab h 


sin 


^m7tx^ 

a J 


sin 


^n7iy^ 


V D J 


The angular natural frequencies are 


<D r 


- 

’ph 




+ 




2\ 


a J 


v D / 


(9.9) 


(9.10) 
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9.1.1.3 Rectangular Plate Example 

Consider a rectangular circuit board made from G10 material. The modulus of elasticity is 2.7e+06 
lbf/in A 2 with Poisson ratio of 0.12. The dimensions of the circuit board are 4 in x 2 in x 0.063 in. The 
board is simply-supported about its perimeter. Assume that the electronic components do not add any 
stiffness. The board has a uniform mass distribution. The total mass is 0.113 Ibm. This includes the G10 
board and the electronic components. The mass per area is 

ph = 0.0141 lbm/in 2 = 3.65e-05 lbf sec 2 /in 3 (9.1 |) 

The circuit board stiffness factor is 


D = 


(2.7e + 06 lbf / in 2 )(0.063 in) 3 

12(l-0.12 2 j 


= 57.1 lbf in 


(9.12) 


The fundamental angular natural frequency is 


®11 = . 


57.1 lbf in 


3.65e-05 lbf sec 2 /in 3 


r \ 


v 4in y 


+ 




v 2in y 


= 3875 rad/sec 


(9.13) 


The fundamental frequency is equivalent to 613.7 Hz. 


The mode shape is shown in the following figure. 


Mode 1 fn= 613.7 Hz 



Figure 9.2. Circuit Board, Fundamental Bending Mode Shape 
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9.1.2 Circular Plate, Bending Vibration 


9.1.2.1 Circular Plate Equations 



Figure 9.3. Circular Plate 


The governing fourth-order spatial equation of motion for the out-of-plane displacement is 


V 4 Z(r,0)-p 4 Z(r,0) = O 


(9.14) 


The biharmonic operator is 


v 4 = V 2 V 2 


d 2 Id 

-^—I-p 

v Sr 2 r Sr 


i s z 

r 2 S6 2 


' d 2 Id 

-r—I-h 

v Sr 2 r Sr 


1 S 

r 2 S0 2 


(9.15) 


The wave number is 


P = 


co 2 ph 

D 


(9.16) 


The kinetic energy in the plate is 


K E = 


phoo 2 r27i fR 2 

2 Jo Jo 


rdrdG 


(9.17) 
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The potential or strain energy in the plate is 


P E 


D r27r j-R 

TJo Jo 


' d 2 Z 1 dZ 

- 7T "I-1" 

v dr 2 r dr 


1 d 2 Z^ 
r 2 d0 2 


2(1-n) 


d 2 Z 

d 2 r 


1 dZ 


-h 

r dr 


1 d 2 z' 

r 2 d0 2 , 


+ 2(1-jli). 


dr 


1 dZ 

7d0 


r drd0 


The solution of equation (9.14) proceeds by assuming separation of the spatial variables. 

Z(r,0) = R(r)©(0) 

Further information on this solution method is given in Reference [8]. 

9.1.2.2 Circular Plate, Simply-Supported 

The boundary conditions for the simply-supported plate are 

Z(a,0) = 0 

M r = 0 at r = a 

The mass-normalized, fundamental mode shape is 

Z 00 = !2^{ J 0 (2.2215 r/a)-0.03686 I 0 (2.2215 r/a)} 


(9.18) 


(9.19) 


(9.20) 

(9.21) 


(9.22) 


Note that I 0 and J 0 in equation (9.22) are Bessel functions. This fundamental mode shape is 
sometimes referred to as an “oil can” mode. The corresponding fundamental frequency is 


®00 - 


4.979 nr 

a 2 VP h 


for p = 0.33 


(9.23) 


A characteristic equation for variable mode numbers and Possion ratios is given in Reference [8]. The 
roots of this equation determine the numerical coefficient on the righthand side of equation (9.23). 
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9.1.2.3 Circular Plate Example, Simply-Supported 

Assume a 48 inch diameter, 0.5 inch thick, aluminum circular plate, with a simply-supported 
circumference. Calculate the fundamental frequency and mode shape. 

The plate stiffness factor is 


(l.0e + 07 lbf /in 2 )(0.5 in) 3 

D = A- - - L -= l. 169e+05 lbf in 

12M-0.33 2 ) 


(9.24) 


The mass per area is 

ph = (0.000259 lbf sec 2 /in 4 )(0.5 in)=1.295e-04 lbf sec 2 /in 3 
The angular natural frequency is 


®00 


4.979 1.169e+05 lbf in 

(24 in) 2 V 1.295e-04 lbf sec 2 /in 3 


= 259.7 rad/sec 


The angular frequency is equivalent to 41.33 Hz. 


(9.25) 


(9.26) 


Model fn= 41.33 Hz 



Figure 9.4. Circular Plate Fundamental Mode Shape 
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9.1.3 Honeycomb Sandwich Plate 

9.1.3.1 Honeycomb Sandwich Plate Characteristics 




Figure 9.5. Sample Floneycomb Sandwich Plates 

ti 



Figure 9.6. Floneycomb Sandwich Geometry 

Honeycomb sandwich structures are designed to have a high stiffness-to-mass ratio. The stiff, strong 
face sheets carry the bending loads, while the core resists shear loads. The face sheets are typically 
made from aluminum or carbon fiber with epoxy resin. The honeycomb core material is usually 
aluminum for aerospace applications. Other core materials include Nomex aramid or Kevlar para- 
aramid fiber sheets saturated with a phenolic resin. In addition, closed cell foams such as Rohacell are 
substituted for honeycomb in some sandwich panel designs. 


Table 9.1. Floneycomb Sandwich Plate Normal Modes 


Domain 

Characteristic 

Low Frequencies 

Bending of the entire structure as if were a thick plate 

Mid Frequencies 

Transverse shear strain in the honeycomb core governs the behavior 

High Frequencies 

The structural skins act in bending as if disconnected 
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The modal response of honeycomb sandwich plates has the characteristic shown in Table 9.1. The 
lower natural frequencies may be calculated using the plate formulas in Sections 9.1.1 and 9.1.2 by 
properly accounting for the material and geometric properties. 

Assume 

1. The skin elastic modulus is much greater than the core modulus. 

2. Each skin has the same material. 

3. Each skin is "thin" relative to the core. 

The honeycomb sandwich plate bending stiffness from Reference [9] is 


D = 


1 

M 

_1 


tj t 2 h 2 



*1 +t 2 


, h-h c + ~[h +t 2 ] 


(9.27) 


The elastic modulus of each skin or face sheet is E f . 


9.1.3.2 Honeycomb Sandwich Plate Example, Part I 

A suborbital rocket vehicle has a circular bulkhead made from honeycomb. The material is aluminum 
for both core and skins. The bulkhead properties are shown in the following table. Calculate the 
fundamental frequency of the bare bulkhead. 


Table 9.2 Bulkhead Parameters 


Parameter 

Value 

Boundary Condition 

Simply Supported 

Diameter 

40 inch 

Core Thickness 

1.0 inch 

Thickness of Each Skin 

0.063 inch 

Total Thickness 

1.125 inch 

Skin Elasticity 

10.0e+06 lbf/in A 2 

Core Elasticity 

Negligible 

Skin Density 

0.10 lbm/in A 3 

Core Density 

0.01 lbm/in A 3 

Poisson's Ratio 

0.3 
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The plate stiffness is calculated in two steps. 


h = 1.0 in + 1[ 0.063 + 0.063] in = 1.063 in 


(9.28) 


D = 


The total skin mass is 


10e + 061bf /in 2 


(0.063in)(0.063in) (1.063in) 2 

1-0.33 2 


0.063in + 0.063in 


= 3.99e + 05 lbf in (9-29) 


m s kin = 2f-l(40in) 2 (0.063 in)(o.l lbm / in 3 ) = 15.83 lbm 


(9.30) 


The core mass is 


ITL. 


f 

n 


v^y 


(40in) 2 (1 in)(0.01 lbm/in 3 j = 12.57 lbm 


The total bulkhead mass is 


m total = 15.83 lbm +12.57 lbm = 28.4 lbm 

The total mass per area is 

28 4 lbm o /-) o 

p A =---— = 0.0226 lbm / in 2 = 5.854e - 05 lbf sec 2 / in 3 

(jt/4)(40in) 2 

The angular fundamental frequency of the bare bulkhead for p=0.33 is 


(9.31) 


(9.32) 


(9.33) 


®00 


4.979 
(20 in) 2 


3.99e + 05 lbf in 
5.854e -05 lbf sec 2 /in 3 


:1028 rad/sec 


(9.34) 


The angular frequency is equivalent to 164 Hz. The mode shape is qualitatively similar to that in Figure 
9.4. 

9.1.3.3 Honeycomb Sandwich Plate Example, Part II 

Now assume that 100 lbm of avionics components are added to the bulkhead and that this nonstructural 
mass is spread uniformly across the bulkhead. Neglect the geometry and stiffness of the avionics. 
Repeat the natural frequency calculation. 
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The total mass per area is 


(100 + 28.4 )lbm , , , 

--= 0.1022 lbm / in 2 = 2.647e - 04 lbf sec 2 / in 3 

(7t/4)(40in) 


The angular fundamental frequency of the mass-loaded bulkhead for |J=0.33 is 


(9.35) 


®00 - 


4.979 
(20 in) 2 


3.99e + 05 lbf in 
2.647e - 04 lbf sec 2 /in 3 


= 483.3 rad/sec 


(9.36) 


The equivalent frequency is 76.9 Hz. 
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10 Damping 

10.1 Damping Metrics 

The following metrics are demonstrated by example throughout this document. The viscous damping 
ratio is defined in terms of the damping coefficients as 

^ = C / C r (10.1) 

The critical damping coefficient is 

c r =2^/km ( 10 . 2 ) 

The amplification factor is 

Q = l'(2y (10.3) 

The amplification factor can be calculated from measured frequency response function data via the half¬ 
power method. 

Q = f /A f (104) 

where f is the peak center frequency and Af is the difference between the two -3 dB points on either 
side of the peak. 

The loss factor is 

11 = 2^ (10.5) 

The logarithmic decrement is 

8 = 2ns/ > /l^ (10.6) 
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10.3 Four Damping Categories 

The four damping types are summarized in the following table. The description is for the free vibration 
of a single-degree-of-freedom system due to initial displacement or velocity. 


Table 10.1. Damping Types 


Type 

Value 

Response Description 

Undamped 

4 = 0 

Simple harmonic motion, sinusoidal response 

Underdamped 

0<^<1 

Damped sine response with exponential decay 

Critically damped 

4=1 

Border between the overdamped and underdamped cases 

Overdamped 

4 >i 

Sum of two decaying exponentials with no oscillation 


The first three of these types are shown by examples in Figure 10.1. The underdamped type 
emphasized throughout this document. 
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Figure 10.1. SDOF Response to Initial Displacement for Three Damping Cases 
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10.4 Damping Mechanisms 

Damping occurs as vibration energy is convert to heat, sound or some other loss mechanism. Damping 
is needed to limit the structural resonant response. Common sources are: 

1. viscous effects 

2. Coulomb damping, dry friction 

3. aerodynamic drag 

4. acoustic radiation 

5. air pumping at joints 

6. boundary damping 

The dominant source for assembled structures is usually joint friction. This damping mechanism may be 
nonlinear due to joint microslip effects. The damping value tends to increase at higher excitation levels. 

Damping may also decrease as the natural frequency increases, such that the amplification factor 
increases with natural frequency. The equation per Reference [10] for approximating Q for an 
electronic system subjected to a sine base input is 


Q = A 


\ 0.76 




0.6 


(10.7) 


where 


A 

= 1.0 for beam-type structures 
= 0.5 for plug-in PCBs or perimeter supported PCBs 
= 0.25 for small electronic chassis or electronic boxes 

Gin 

Sine Base Input (G) 


Beam structures: several electronic components with some interconnecting wires or cables. 

PCB: printed circuit board well-populated with an assortment of electronic components. 

Small electronic chassis: 8-30 inches in its longest dimension, with a bolted cover to provide access 
to various types of electronic components such as PCBs, harnesses, cables, and connectors. 
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10.5 Huntsville, Alabama Pedestrian Bridge Damping 



Figure 10.2. University Drive Pedestrian Bridge 



Figure 10.3. Slam Stick X, Triaxial Accelerometer & Data Logger, Shown on Book Shelf 

The bridge in Figure 10.2 is near the University of Alabama in Huntsville. The author walked to the 
center of the bridge and mounted the triaxial accelerometer in Figure 10.3 on the deck floor using 
double-sided tape. He was the only person on the bridge. The ambient vibration response of the bridge 
was negligibly low in each of the axes, while he remained standing. The wind was very light on this day. 
The author then jumped up-and-down on the bridge to excite its vibration modes. The response in the 
vertical axis was significant. The response in each of the lateral axes remained negligibly low. 
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Figure 10.4. Signal Identification via Damped Sine Curve-fit 

The jumps were performed near the 80 second mark. The accelerometer data was band-passed filtered 
from 1 to 20 Hz. The response after 82 seconds represents the free vibration decay. A damped sine 
curve-fit synthesis was performed on the vertical acceleration time history as shown in Figure 10.4, 
using trial-and-error with convergence. The bridge’s fundamental frequency is 2.2 Hz with 0.16% 
damping. The damping ratio is very low. But note that the fundamental frequency and damping may be 
nonlinear. A pedestrian’s vertical forces correspond to each footfall, and typically occur at 2.0 Hz. This 
is very close to the bridge’s 2.2 Hz natural frequency. 

Also note that there is a potential for pedestrians to synchronizing their steps with the bridge motion 
and with one another. This behavior is instinctive rather than deliberate. Pedestrians find that walking 
in synchronization with the motion of a bridge is more comfortable, even if the oscillation amplitude is 
initially very small. This cadence makes their interaction with the movement of the bridge more 
predictable and helps them maintain their balance. But the synchronization also causes the pedestrians’ 
gait to reinforce the bridge’s oscillation in a resonance-like manner. These sorts of problems occurred 
after the opening of the London Millennium Bridge in 2000. Both passive and tuned mass dampers were 
added to the bridge for vibration control. Retrofitting the University Drive Bridge with dampers is 
unnecessary due to its low pedestrian traffic volume. But this could be an interesting project for the 
nearby engineering students. 


Acceleration Time History University Drive Bridge 

“i-1-1-1-1-1-1- r~ 


Measured Data 
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10.6 Pegasus Launch Vehicle Damping 

10.6.1 Drop Transient 



Figure 10.5. Pegasus Launch Vehicle, Drop & Stage 1 Burn 



Figure 10.6. Pegasus Fundamental Bending Mode Shape, Exaggerated 

A modified L-1011 aircraft carries the Pegasus vehicle up to an altitude of nearly 40,000 ft and a speed of 
Mach 0.8, as shown in Figure 10.5. Pegasus is suspended underneath the aircraft by hooks, where it 
develops an initial displacement due to gravity. The strain energy is suddenly released at the onset of 
the drop transient, causing Pegasus to oscillate nearly as a free-free beam. This is a significant “coupled 
loads” event for the payload which is enclosed in the fairing at the front end of the vehicle. The Pegasus 
first stage then ignites. The payload is eventually delivered into a low earth orbit. A certain 
Pegasus/payload configuration was analyzed via a finite element modal analysis. The resulting wire mesh 
model of the fundamental bending mode is shown in Figure 10.6. Note that this was a different Pegasus 
configuration than that represented by the flight data in Figure 10.7. 
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10.6.2 Drop Transient Damped Sine Curve-fit 



Figure 10.7. Pegasus Drop Transient, Flight Accelerometer Data, Free Vibration Response 

The data was measured in the transverse axis at the payload interface. The response is nearly a 
textbook quality damped sinusoid, with an exponential decay. The natural frequency and damping ratio 
are identified via a synthesized damped sine, curve-fit. 
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10.6.3 Drop Transient Logarithmic Decrement Method 



Figure 10.8. Pegasus Drop Transient, Flight Accelerometer Data, Logarithmic Decrement 
The logarithmic decrement for a starting peak and a peak n cycles later is 


8 = 




n 




( 10 . 8 ) 


The logarithmic decrement for the two peaks shown in Figure 10.8 is 


8 = 


If 0.7473 
8 n { 0.4026, 


0.0773 


(10.9) 


The logarithmic decrement value is equivalent to 1.2% damping in agreement with the previous damped 
sine curve-fit. The logarithmic decrement method has been included for historical reasons. The 
damped sine curve-fit method is more robust. 
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10.7 Transamerica Building Damping 



Figure 10.9. Transamerica Pyramid 


The Transamerica Pyramid is built from a steel frame, with a truss system at the base. The height is 850 
ft (260 m). Reference [II] gives natural frequency and damping as obtained in the 1989 Loma Prieta 
earthquake and from ambient vibration. The ambient vibration was presumably due to wind, low level 
micro-tremors, mechanical equipment, outside street traffic, etc. 


Table 10.2. Transamerica Pyramid, Modal Parameters 


Direction 

Loma Prieta Earthquake 

Ambient Vibration 

fn (Hz) 

Damping 

fn (Hz) 

Damping 

North-South 

0.28 

4.9% 

0.34 

0.8% 

East-West 

0.28 

2.2% 

0.32 

1.4% 


The results show non-linear behavior with an increase in damping during the severe earthquake relative 
to the benign ambient vibration. 
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11 Sinusoidal Vibration 

11.1 Descriptive Statistics 

SENE FUNCTION 



Figure 11.1. Sine Function, One-Cycle 

A sine function is characterized by its frequency and amplitude. The frequency is the inverse of the 
period, where the period is the duration of one cycle. The amplitude may be represented by any of the 
three parameters shown on the right-hand-side. The zero-to-peak and RMS amplitudes for sinusoidal 
oscillations are related by 

A Peak = Arms (11.1) 

Note that 

Arms 2 = A mea n 2 + A lcj 2 (I 1.2) 

where is the standard deviation 
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The RMS value is thus equation to the standard deviation for the case of zero mean. A pure sine 
function always has zero mean. 


Time History 
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Counts 


Figure 11.2. Sine Function & Flistogram 

A sine function and its corresponding histogram are shown in Figure 11.2. The histogram divides the 
amplitude points into bins and then counts the points in each bin. The number of points in each bin are 
then plotted as a bar graph. Some of the histograms in this document are turned 90 degrees clockwise 
and shown adjacent to the corresponding time history for educational purposes. 

A sine function has a bathtub histogram. The amplitude tends to dwell at either the peak or valley with 
a fairly rapid transition in between. Histograms are important for fatigue analysis. 

Now let x(t) be a sinusoidal displacement function. The displacement equation is 


x(t) = Xsin(cot) 


(11.3) 


where 


X is the zero-to-peak amplitude 
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The corresponding velocity v(t) is 
v(t) = ooXcos(oot) 

The acceleration a(t) is 

a(t) = -(o 2 Xsin((ot) 

Table 11.1. Peak Values Referenced to Peak Displacement 


Parameter 

Value 

displacement 

X 

velocity 

oo X 

acceleration 

oo 2 X 


Table 11.2. Peak Values Referenced to Peak Acceleration 


Parameter 

Value 

acceleration 

A 

velocity 

A / co 

displacement 

A / co 2 


(11.4) 


(11.5) 


Let A be the peak acceleration. 

A = co 2 X (11.6) 

Conversely, 


X = A / co 2 (11.7) 

11.2 Sine Testing 

Table I 1.2 has enormous consequences for vibration testing. Commercial shaker tables, such as those 
found in Figure 4.9 and Figure 4.10, are typically limited to a few inches peak-to-peak displacement 
unless they are large-scale tables specifically designed for earthquake simulation. Vibration tests are 
usually specified in terms of acceleration. The test specification starting frequency is typically 10 to 20 
Hz to limit the displacement. 
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11.3 Sample Sine Calculations 

Astronaut Michael Collins wrote [12]: 

The first stage of the Titan II vibrated longitudinally, so that someone riding on it would be 
bounced up and down as if on a pogo stick. The vibration was at a relatively high frequency, 
about 11 cycles per second, with an amplitude of plus or minus 5 Gs in the worst case. 

The equivalent angular frequency for II Hz is 69.1 rad/sec. Calculate the corresponding displacement X 
using equation (I 1.7). 

X = (5 G )(386 in/sec 2 /G)/( 69.1 rad/sec) 2 = 0.40 inch zero-to-peak (11-8) 

M.J. Griffin [13] cited a study where 5% of the passengers on a ship experienced vomiting while the ship 
underwent a 0.2 m/sec A 2 acceleration amplitude during a 3 hour. Assume a 0.2 Hz frequency (1.26 
rad/sec A 2). Calculate the displacement. 

X = (0.2 m/sec^2)(100 cm/m) / ( 1.26 rad/sec) 2 = 12.7 cm zero-to-peak (11.9) 

11.4 Beat Frequency Effect 


Channel Beam Free Vibration 



Figure 11.3. Channel Beam Free Vibration Response at Free End to Force Impulse 

Consider the case of two closely-space vibration frequencies generating sound. The envelope of the 
two signals oscillates at a frequency which is the differences of the two frequencies. This beat frequency 
is detected by the human hearing and brain system, but it is a psychoacoustic effect which does not 
affect the mechanical hardware. 
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A channel beam was mounted in a cantilever configuration. The free end was struck with an impulse 
hammer. The resulting acceleration time history at the free end is shown in Figure 11.3. A synthesized 
damped sine curve-fit yielded natural frequencies at 583 and 691 Hz. The beat frequency is the 
difference which is 108 Hz. 



Figure 11.4. Channel Beam Free Vibration Response Fourier Transform 

A Fourier transform of the beat frequency data is shown in Figure 11.4. Peaks occurs 583 and 691 Hz 
but not at the 108 Hz beat frequency. Fourier transforms are covered in Section 15. 
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12 Dynamic Models for External Excitation 

12.1 Two models 

There are two basic SDOF spring-mass models used for external excitation in shock and vibration 
analysis. One is the applied force model, and the other is the base excitation model. An introduction to 
these types and images was previously given. The purpose of this section is to demonstrate solution 
methods for each of these type in both the time and frequency domains. 

12.2 Applied Force 

12.2.1 SDOF Response Equations 

The applied force model and its free-body diagram are 



f 



x 


Figure 12.1. Single-degree-of-freedom Model, Applied Force 
The equation of motion for applied force F is derived using Newton’s law. 

p = mx 

mx = f(t)-kx-cx 

m x + cx + kx = f (t) 

2 1 

x + 2£,co n x + co n x = —f(t) 
m 


( 12 . 1 ) 

( 12 . 2 ) 


(12.3) 

(12.4) 
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The resulting equation of motion is non-homogenous. It can be solved via Laplace transforms if the 
force is deterministic such as a sine function. A convolution integral is needed if the force varies 
arbitrary with time. The resulting displacement for an applied force is 


x(t)= —f(x) |exp[-^co n (t—x)]|{sin(o d (t-x)}ch 


mco d 


(12.5) 


The convolution integral in equation (12.5) can be solved in terms of a ramp invariant digital recursive 
relationship [14]. 

Alternate methods for numerically solving equation (12.4) include the Runge-Kutta and Newmark-beta 
methods. The Runge-Kutta method may unstable for stiff systems, however. 

The displacement convolution integral has an embedded impulse response function. 


h d (t)= —[exp(-£ra n t)][sm(a d t)] 


mco d 


( 12 . 6 ) 


The Laplace transform of the impulse response function is the steady-state frequency response function 
(FRF). 


H d (s) = — 
m 


1 

2 2 

s +2^co n s + co n 


(12.7) 


This FRF can be represented in the frequency domain by setting s = jco, where 


j = \[~i 


H d (®) = 


l 

m 



( 12 . 8 ) 


The FRF in equation (12.8) is the Fourier transform of the impulse response function in equation (12.7). 
Fourier transforms are covered in Section 15. An equivalent FRF is 


H d (a>) = 


k 




(12.9) 
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12.2.3 Input & Response Example, Applied Force 

The single-degree-of-freedom system in Figure 12.1 has a mass of 20 kg, stiffness of 78960 N/m, natural 
frequency of 10 Hz and an amplification factor Q=10 for this example. It is excited into resonance by a 
20 N applied force sinusoidal excitation at 10 Hz. The exact response can be calculated via a Laplace 
transform solution. The result is shown in Figure 12.2. 




Figure 12.2. Sine Applied Force & Response 


The peak response is 2.53 mm. Note that the following ratio is approximately equal to the Q value. 

( 12 . 10 ) 


kx 


(78960 N/m)(0.00253 m) 

F 


20 N 


10 
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0.1 1 10 



Figure 12.3. Transfer Function for Steady-State Force Excitation 

The response calculation could be repeated for numerous frequency and amplification factor cases to 
determine the maximum steady-state response for each case. The transfer function represents the 
displacement divided by the applied force as normalized by the stiffness such that the magnitude is 
unitless, as previously introduced. The transfer function magnitude is 


kx 

F 


_I_ 

^(ivf+^p) 2 


p=f/f n 


( 12 . 11 ) 


The phase angle by which the displacement lags the force is 


9 = arc tan 


2^p/(l-p 2 ) 


( 12 . 12 ) 


The applied force transfer function in Figure 12.3 can roughly be divided into three zones. The function 
is controlled by stiffness at the low frequency ratio end. Its magnitude converges to the simple static 
Hooke’s law as the frequency ratio approaches zero. The function is controlled by damping for 
frequency ratios at or near the resonance condition where the ratio is equal to 1. The function is 
controlled by mass at higher ratios, above V2. 
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12.2.4 Six Transfer Functions for Applied Force and Response 

Recall equation (12.8). The following equation represents the ratio of displacement to force in terms of 
their respective Fourier transforms. 


x(<») 

F(«0 


H d (o>) 


( 12 . 13 ) 


Similar frequency response ratios can be formed for the velocity and acceleration responses relative to 
force. Inverses of these functions can also be taken. The names of these transfer functions are given in 
the following tables. 


Table 12.1. Transfer Function Nomenclature 


Displacement / Force 

Velocity / Force 

Acceleration / Force 

Admittance, 

Compliance, 

Receptance 

Mobility 

Accelerance, 

Inertance 


Table 12.2. Reciprocal Transfer Function Nomenclature 


Force / Displacement 

Force / Velocity 

Force / Acceleration 

Dynamic Stiffness 

Mechanical Impedance 

Apparent Mass, 

Dynamic Mass 


Furthermore, the receptance frequency response function can be modeled in terms of a block diagram. 
Similar diagrams could be drawn for the other transfer functions. 


Input Force 

F H 


> 


Transfer Function 
Hd(<») 


Displacement Response 

-► 

X(to) 


Figure 12.4. FRF Diagram with Input & Output 
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12.2.5 Helicopter Ground Resonance 




Figure 12.5. Helicopter Preparing for Takeoff and SDOF System Analogy 


A sketch of a generic helicopter is shown in Figure 12.5. The helicopter may be modeled as a single- 
degree-of-freedom system excited by an applied force. 

A new helicopter design must undergo thorough testing to ensure the reliability of the design with 
respect to vibration. It may encounter severe vibration while it is on the ground, preparing for takeoff. 
A similar problem may occur immediately after landing. 

Many helicopters have fully-articulated rotor systems. The advancing blades flap upward and the 
retreating blades flap downward in order to maintain uniform lift during forward flight. A helicopter with 
a fully-articulated rotor has some measure of an unbalanced rotational force, which causes an oscillating 
force. This force may excite one of the helicopter’s natural frequencies, as the rotor accelerates to its 
full operating speed prior to takeoff. This condition is called resonant excitation. 

John O. Emmerson wrote in Reference [15]: 

In the mid-1950s, the emerging helicopter introduced a frightening problem called ground 
resonance. A perfectly sound production helicopter firmly supported on its landing gear with 
the rotor turning at or near flight rpm would suddenly initiate a “tramping” on the gear and 
then, almost immediately, completely disintegrate. The Sikorsky company, which was leading 
the pack in developing the helicopter, quickly gathered the most data on these ground-based 
disasters and thus shouldered the lion’s share of corrective understanding and action. The 
company established that the lead/lag pivots on the rotor blades, the pylon elastic stiffness, and 
the elastic components of the landing gear, including the tire and shock struts, all working 
together as a single structural entity, had an unfortunate resonant frequency near the rotor 
operating rpm. The engineers were able to show that a fearless pilot, entering the resonant 
phase, could save the day be pulling up the collective pitch control and unloading the landing 
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gear (i.e., getting airborne). Lowering the resonant frequency and adding damping offered a 
permanent solution. 

Today, all of us in the helicopter development business, with the apprehension borne of 
substantial ignorance, approach initial ground whirl testing of each new helicopter with 
tremulous respect for the demon of ground resonances. 


Modern helicopters avoid ground resonance by using dampers on the blades and on the landing gear. 
The landing gear may have shock absorbing struts. 



Figure 12.7. TH-55 Osage, Destroyed by Ground Resonance 

An intact TH-55 helicopter is shown in Figure 12.6. A similar TH-55 helicopter was destroyed by 
ground resonance in March 1967, at Ft. Wolters, TX, as shown in Figure 12.7. 
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12.2.6 Two-Degree-of-Freedom System Subjected to Applied Forces 
12.2.6.1 Two-Degree-of-Freedom System Equations of Motion, Time Domain 



x 2 



m, 

_| x, - k 2( x |- x 2> 

■ f 

J 

1 

L 




T 

/ / / / / 



c 2 (ii -x 2 ) 
x i 


_k l x l -C|X I 


Figure 12.8. Two-DOF System with External Forces 


Recall the grounded, two-DOF system from Section 6.1.1. It is shown again in the above figure with 
added forces. The coupled equations of motion are 


m i 

o “ 

"V 

+ 

C 1 + c 2 

~ c 2~ 

’ x l’ 

+ 

'ki+k 2 

- k 2 _ 

’ x l’ 

_ 


0 

m 2_ 

_ x 2_ 


_ —c 2 

c 2 . 

_ x 2_ 


~ k 2 

k 2 . 

_ x 2_ 


_ f 2_ 


Note that each of the coefficient matrices is symmetric. This is typical for linear structural dynamics 
systems. 
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A shorthand form of equation (12.4) is 


Mx+Cx+Kx=f 

The matrices and vectors are 


(12.15) 


M = 

m l 

0 ~ 

, c = 

c l +c 2 ~ c 2 


0 

m 2_ 


_ “ c 2 c 2 _ 


ki +k 2 

<N 

1 

_ 

fl 


X 1 

—k 2 

k 2 . 

, f = 

f 2 

, X = 

x 2 


(12.16) 


The natural frequencies and mode shapes are found from the undamped, homogeneous equation as was 
shown in Section 6.1.1. Now define a modal coordinate T|(t) in terms of the normalized eigenvector 

matrix Q such that the displacement vector is 


x(t) = Qn(t) 


~qn 

A21 


4l2 

<122 


^(t) 


Repeat the steps in Section 6.1.1 with the added force term. 


Q T MQii + Q T CQii + Q T KQrj = Q T f 


(12.17) 


(12.18) 


The orthogonality relationships yield 


I f) + Q t CQi) + Qr[ = Q t F 


(12.19) 


The equation of motion can now be written with the modal damping simplification as 


1 O' 

V 

-f 

’2^1©! 

0 

m’ 

+ 

COj 2 

0 

"ni’ 


qn 

q2i 

'fi 

0 1 

J\2. 


0 

2^2 co 2_ 

.02_ 


0 

-1 

<N 

<N 

3 

_ r l2_ 


_qi2 

q22. 

_ f 2 


( 12 . 20 ) 


The two equations are now uncoupled in terms of the modal coordinates. Now consider the case 
where the forcing functions are deterministic time domain functions. The modal displacement for the 
response to the harmonic forces can be found via Laplace transforms. Alas, the resulting modal 
displacement equations even for sinusoidal forces are too unwieldy for inclusion in this book but are 
given in Reference [16]. The physical displacements are then found via equation (12.17). The solution 
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for the case of arbitrary forces can be calculated using the ramp invariant digital recursive filtering 
relationship in Reference [14]. 

12.2.6.2 Two-Degree-of-Freedom System, Time Domain Example 

The system in Figure 12.18 has the parameters previously shown in Table 6.1, with 5% damping for 
each mode. Its natural frequencies are again 20.0 and 42.5 Hz. The mass-normalized mode shapes in 
matrix format are 


8.93 - 10.64 

15.05 12.63 


( 12 . 21 ) 


It is subjected to a single sinusoidal force applied to the top mass, with zero initial conditions. The 
amplitude is 1 Ibf. The forcing frequency is 20.0 Hz, which is the same as the fundamental frequency. 
The results are shown in Figure 12.9. 

Now consider the steady-state amplitude ratios between the mass responses. The responses are very 
nearly in proportion to the eigenvector coefficients in the first column of equation (12.21). The 
response of mass 2 is thus nearly (15.05/8.93) times greater than that of mass 1. The measured ratio is 
slightly inexact because the 20.0 Hz excitation also excites the 42.5 Hz mode, but the second mode’s 
participation can be considered as negligibly low due to the greater than one-octave frequency 
separation. 

The analysis is repeated with the forcing frequency at 42.5 Hz, as shown in Figure 12.10. The steady- 
state response of mass 2 is nearly (12.63/10.64) times than that of mass 1, per the second column of 
equation (12.21). The negative sign in the second column indicates that the responses are 180 degrees 
out-of-phase with one another. 
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Figure 12.9. Displacement Responses to 20.0 Hz Excitation 



Figure 12.10. Displacement Responses to 42.5 Hz Excitation 
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12.2.6.3 Two-Degree-of-Freedom System Equations of Motion, Frequency Domain 

Equation (12.20) can be used in the frequency domain for the case of steady-state harmonic forces. 
Define force and modal displacement coordinates as follows 


f (t) = Fexp(jcot) (12.22) 

r[(t) = Nexp(jcot) (12.23) 


Substitute these harmonic terms into Equation (12.20). This process is akin to taking a Fourier 
transform of the equation. The result for modal displacement i after simplification is 


co; -co 


+ j2^coico 


Ni=qiiFi+q2iF 2 


(12.24) 


Nj = 


qii F i+q 2 i F 2 



+ j2^ i C0 i G) 


(12.25) 


The modal displacement equation can be generalized to a system with n degrees-of-freedom as 


Ni = I 

p=l 


q P iFp 


coj 2 -co 2 ) + j2£, j cojco 


The physical displacements are then found from 


X(cd) = QN(g)) = 


qn qi2 
q2i q22 


N(co) 


The physical displacement at degree-of-freedom k is 


(12.26) 


(12.27) 


x k = Iqti N i ( |2 - 28 ) 

i=l 

The principle of linear superposition holds in this case. The displacement can be calculated for an 
individual force. Then the total displacement is the sum of the displacements for the individual forces. 
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Assume that the only force acting on the system is at degree-of-freedom p The modal displacement at 
degree-of-freedom i is 


Ni = 


^pi^p 


coj 2 -oo 2 ) + j2£,jcoj(o 


The physical displacement at degree-of-freedom k for the single force case is 


( 12 . 29 ) 


x k =I 

i=l 


qkiq p iFp 


^(0j 2 -to 2 j +j2^j(0j(0 

The receptance function for the single force applied becomes 


( 12 . 30 ) 





E< 

i=l 


^kiQpi 


CDj 


co + j2£ ) j coj co 


( 12 . 31 ) 


The numerator in the series shows the relationship between the receptance function and the mass- 
normalized eigenvector terms. 
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12.2.7 Two-Degree-of-Freedom System, Frequency Domain Example 



Frequency (Hz) 

Figure 12.11. Two-DOF System Receptance FRF 

The system in has the parameters previously shown in Table 6.1, with 5% damping for each mode. Its 
natural frequencies are again 20.0 and 42.5 Hz. Its receptance frequency response function from 
equation (12.31) is given in Figure 12.11 for four permutations. The first index in the legend is the 
response location. The second is the applied force location, where each location is a degree-of- 
freedom. The phase angles could also be plotted but are omitted for brevity. 
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12.3 Base Excitation 

12.3.1 SDOF Response Equations 

This case is also referred to as support motion or seismic excitation. The base excitation model and its 
free-body diagram are 




k(y-x) c(y-x) 


Figure 12.12. Single-degree-of-freedom Model, Base Excitation 

The base displacement is y, and base velocity is y. The equation of motion for base excitation is 
derived using Newton’s law. 


XF = mx ( I2.32) 

mx = -k(y-x)-c(y-x) (12.33) 


mx + cx + kx = cy + ky 


(12.34) 


The challenge with equation (12.34) is that the base displacement and base velocity are usually unknown, 
but the base acceleration is known. The workaround is to define a relative displacement z as 


z = x-y 


(12.35) 


Through substitution and algebraic manipulation, the following ordinary differential equation is derived in 
terms of the relative displacement. 


mz + cz + kz = -my 


(12.36) 


The resulting equation of motion is non-homogenous. The my term on the right-hand side can be 
regarded as a virtual force. 
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The equation for base excitation can be further simplified as 

z + 2^G3 n z + co n 2 z = — y (12.37) 

The absolute acceleration is then calculated as 

* = ^ + y (12.38) 

Equation (12.37) can be solved via Laplace transforms if the base acceleration is deterministic such as a 
sine function. A convolution integral is needed if the excitation varies arbitrarily with time. The 
resulting relative displacement equation for acceleration base excitation is 


z(t)= —y(x) {exp[-^co n (t-x)]}{sinco d (t-x)}dT 


(12.39) 


Equation (12.37) can also be solved in the frequency domain. Time and frequency domain analyses 
should yield the same results per Parseval's theorem. 

The relative displacement convolution integral has an embedded impulse response function. 

h z (t)=—[exp(-£© n t)] [sin (eo d t)] (12.40) 

The Laplace transform of the impulse response function is the steady-state transmissibility function for 
relative displacement. 


H z (s) = 


1 

2 2 

s + 2^co n s + co n 


(12.41) 


This transmissibility function can be represented in the frequency domain by setting s = jco, where 



H z (co) = 



(12.42) 


The relative displacement transmissibility function in equation (12.42) is the Fourier transform of the 
impulse response function in equation (12.40). Fourier transforms are covered in Section 15. 
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The absolute acceleration response convolution integral is 


x(t) = — |q yCO {exp[-£<*> n (t - x)]J {(co d 2 - £, 2 a> n 2 jsin co d (t - x) + 2^co n co d cos co d (t - x)jdx 


(12.43) 

The convolution integral in equation (12.43) is computationally inefficient, however. An alternative is to 
use the Smallwood ramp invariant digital recursive filtering relationship [17]. The recursive filtering 
algorithm is fast and is the numerical engine used in almost all shock response spectrum software. It is 
also accurate assuming that the data has a sufficiently high sample rate and is free from aliasing. The 
algorithm is discussed further in Section 18.5. 

The embedded impulse response function in equation (12.43) is 


h a(t) = 


J_r 

co d 


exp(-%co n t) {ca d 2 -^ 2 co n 2 jsin(co d t) + 2^(o n co d cos((o d t) 


(12.44) 


The Laplace transform of the impulse response function is the steady-state transmissibility function for 
absolute acceleration. 


H a (s) = 


9 

2^co n s + co n 

2 2 

s +2^co n s + co n 


The transmissibility function in the frequency domain is 


(12.45) 


H a (co) = 


co n z +j2^co n co 
{co n 2 — co 2 ) + j 2^co n (o 


(12.46) 


12.3.2 Input & Response Example, Base Excitation 

The single-degree-of-freedom system in Figure 12.12 has a natural frequency of 10 Hz and an 
amplification factor Q=10 for this example. It is excited into resonance by a 1 G sinusoidal base 
excitation at 10 Hz. The exact response can be calculated via a Laplace transform solution. The result 
is shown in Figure 12.13. 
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Figure 12.13. Sine Base Input & Response 
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Figure 12.14. Transmissibility Function for Steady-State Sine Base Input 

The response calculation could be repeated for numerous frequency and amplification factor cases to 
determine the maximum steady-state response for each case. The transmissibility results for three 
amplification factor cases are shown in Figure 12.14. The frequency ratio is the excitation frequency f 
divided by the natural frequency fn. The magnitude is the steady response divided by the input 
amplitude. The transmissibility tends to unity gain when the excitation frequency is much less than the 
natural frequency. Resonant amplification occurs when the frequency ratio is at or very near one. 
Isolation occurs when the excitation frequency is greater than V2 times the natural frequency. 

Note that the Q value is approximately equal to the peak transmissibility at resonance in Figure 12.14. 
Also, the resonant time domain response in Figure 12.13 reached a steady-state peak of Q times the 
base input. This suggests that sine vibration testing could be a way for measuring the amplification 
factor and its corresponding damping ratio by noting the peak response either in the frequency or time 
domain. But this simple relationship should be considered only for the special case of a single-degree-of- 
freedom system. Multi-degree-of-freedom systems and continuous systems have responses which 
depend on damping as well as on mode shape parameters. 
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Consider a cantilever beam subjected to base excitation exciting its fundamental bending mode. An 
accelerometer is mounted on the beam’s free end, and another at its midpoint. The peak response at 
free end will be higher than that at the midpoint due to the respective location-dependent modal gain 
factors. The half-power bandwidth method is a much more reliable damping identification method for 
measured response data from this beam and from general structures. This method would yield the 
same damping value from either of the two accelerometer data results. See the examples in Sections 
13.4 and 13.5. 

The transmissibility function plotted in Figure 12.14 was introduced in Section 12.3. The magnitude of 
H(p) can be represented as 


| H (P)| = 


1 + (2^P f 


!(l-p 2 ) +(2^p f 


P = f/fn 


(12.47) 


The phase angle by which the response lags the input is 


c|) = arctan 



(12.48) 


12.3.3 Alternate Base Excitation Method using Seismic Mass 


k 



x 2 


X 1 


Figure 12.15. Seismic Mass for Base Excitation, Semi-definite System 

Recall the semi-definite system in Section 6.1.6 and the two-DOF system subjected to an applied force 
in Section 12.2.6. A hybrid of these systems can be configured for base excitation by setting the 
bottom mass to a value vastly larger than the top mass by a factor of, say, le+04. Then an external 
force is applied to the base to produce the desired acceleration. This approach is unnecessary for 
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SDOF system base excitation but becomes very useful for complex, multi-degree-of-freedom systems. 
The equation of motion is 



0 " 

"V 

+ 

c 

-c 

~ x l~ 

+ 

1 

FT 

1 

FT 

_i 

~ x i~ 


"f" 


y" 

_ 0 

m 2_ 

_ x 2_ 


-c 

c 

>2. 


-1 

M 

M 

1 

_i 

_ x 2_ 


_0_ 


L 0 


(12.49) 


The hybrid equation is solved using the methods in Section 12.2.6. There is a rigid-body mode at zero 
frequency. The second modal frequency is 


co 2 = 


i k ( 


in i + nv 


m l m 2 


( 12 . 50 ) 


The second modal frequency becomes 

co 2 ~ yj k/m 2 , for m 1 »m 2 


( 12 . 51 ) 


The second frequency is essentially equal to the top mass and the spring as a single-degree-of-freedom 
system. 

Repeat the example in Section 12.3.2. The natural frequency was 10 Hz in this example with Q=10. 
Use of the seismic mass method requires explicit knowledge of the mass and stiffness values, as shown 
for this example in Table 12.3. The top mass is excited into resonance by a 1 G sinusoidal base 
excitation at 10 Hz via an equivalent 100,000 Ibf force applied to the bottom, seismic mass. The 
acceleration response is shown in Figure 12.16. The result is the same as that in Figure 12.13 with the 
response time histories in both plots reaching 10 G steady-state. 


Table 12.3. Seismic Mass, Hybrid System, Parameters 


Variable 

Value 

mi 

l.0e+05 Ibm 

m 2 

10 Ibm 

k i 

102.2 Ibf/in 
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Figure 12.16. Sine Base Input & Response, Seismic Mass Method 

Note that in some cases it may be necessary to ground the seismic mass with a very soft spring for 
numerical stability purposes, such that the fundamental frequency is just slightly above the zero 
frequency but is still much less than the frequency in equation (12.51). 
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12.3.4 Two-DOF Response Equations, Time Domain 



-c,(x,-y) 


Figure 12.17. Two-DOF System Subjected to Base Excitation 

Consider the system in the above figure, which is sometimes referred to as two-stage isolation per 
Reference [18]. The equation of motion can be derived using Newton’s law. 


ni| 0 

" x l" 


~ C 1 + c 2 

- c 2~ 

"*f 


ki+k 2 

- k 2 _ 

X 1 


c i y+ k i y 

<N 

a 

o 

_ x 2_ 


_ -c 2 

c 2 . 

>2. 


-k 2 

k 2 . 

_ x 2_ 


0 


The familiar challenge with equation (12.52) is that the base displacement and base velocity are usually 
unknown, but the base acceleration is known. The workaround is to define a relative displacement z 
such that 


z i = x i-y (12.53) 

Substitution of the relative displacements into equation (12.52) yields the following equation after 
algebraic manipulation. 
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m l 

0 “ 

z i 


" c l+ c 2 

-c 2~ 

z i 


ki+k 2 

“ k 2 

~ z l~ 


-miy 


0 

m 2 _ 

_ z 2_ 

+ 

. _c 2 

c 2 . 

>2. 

+ 

-k 2 

k 2 . 

_ z 2_ 


_-m 2 y_ 

(12.54) 


The natural frequencies and modes shapes are calculated using the method from Section 6.1.2. The 
modal decoupling process is essentially the same as that in Section 12.2.6, with force vector changed to 


F = 


-miy 

-m 2 y 


(12.55) 


The relative displacement is represented in terms of the mass-normalized eigenvector matrix and the 
modal displacement. 


Z = Q r| = 


’qn 

_§2i 


qi2 

§22. 


o(t) 


The uncoupled equations of motion are 


(12.56) 


1 O' 

"m" 

+ 

"2^i ©i 

0 

2^2 co 2_ 

~v 

+ 

COj 2 

0 

~v 


~qn 

§21 

-miy 

o 

Jl2_ 


0 

Ji2_ 


0 

-1 

(N 

<N 

s 

Jl2_ 


_qi2 

§22 _ 

_-m 2 y_ 


The modal displacement can be calculated via Laplace transforms for a deterministic base acceleration. 
Modal accelerations can be calculated by twice differentiating the modal displacements. 

A ramp invariant digital filtering relationship can be used for a base excitation which varies arbitrarily 
with time. Such equations are available for direct calculation of both the modal displacement and 
acceleration. The filtering equations are taken from Reference [14] given the virtual force and 
displacement characteristics of equation (12.57), rather than the SDOF formulas for base excitation 
given in References [17], [19]. 

The physical relative displacement can then be calculated by applying the modal displacements to 
equation (12.56). The physical relative acceleration can then be determined from the modal 
acceleration via the same eigenvector matrix transformation in equation (12.56). The physical absolute 
acceleration is then found by adding the base acceleration to the physical relative acceleration. 

12.3.5 Two-DOF Response Equations, Time Domain Example 

The system in Figure 12.17 has the parameters previously shown in Table 6.1 with 5% damping for each 
mode. Its natural frequencies are again 20.0 and 42.5 Hz. It is subjected to sinusoidal base excitation, 
with zero initial conditions. The base acceleration amplitude is 1 G with frequency at 20.0 Hz, which is 
the same as the fundamental frequency. The results are shown in Figure 12.18. The analysis is repeated 
with the excitation frequency at 42.5 Hz, as shown in Figure 12.19. 
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Figure 12.18. Acceleration Responses to 20.0 Hz Excitation 



Figure 12.19. Acceleration Responses to 42.5 Hz Excitation 
The responses are 180 degrees out-of-phase with one another. 
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12.3.6 Two-DOF Response Equations, Frequency Domain 

Equation (12.57) can be used in the frequency domain for the case of steady-state harmonic base 
excitation. Define base acceleration and modal displacement coordinates as follows 

y( t ) = YA ex p(j° )t ) 

^(t) = Nexp(jcot) 

Substitute these harmonic terms into Equation (12.57). The result for modal displacement i after 
simplification is 


(12.58) 

(12.59) 


COj -CO ) + j2£,j CGjCG 


N 


*=-( 


qii m i + q 2 i m 2 


) y A 


- (qiimi + q 2 im 2 ) Ya 

Nj = —y-r--- 

I coj 2 -co 2 ) + J2£,j cojco 

The physical relative displacements are then found from 


Z(eo) = Q N(co) 


4ii 

_q2i 


qi2 

q22 


N(co) 


(12.60) 

(12.61) 


(12.62) 


The physical relative displacement at degree-of-freedom k is 


2 2 

Zk = IqkiNi = -Y A £ 


i=l 


i=l 


qki(qn m i + q2i m 2) 


r\ r\ 

coj -co ) +j2^jcojco 


(12.63) 


The transmissibility function relating physical relative displacement to base excitation is 


H r d (®) : 


Zk_ 

y a 


2 

■z 

i=l 


Qki (qii m l + q2i m 2) 


coj 2 -co 2 1 + j2^j cojco 


(12.64) 


The transmissibility function relating physical relative acceleration to base excitation is 


H ra (G>)=co 2 £ 
i=l 




qki qii m i + q2i m 2 


i m 2) 


coj 2 - co 2 ) + j2£,j cojco 


(12.65) 


The transmissibility function relating physical absolute acceleration to base excitation is 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


103 

















enDAQ 


Chapter: Dynamic Models for External Excitation 


Haa (®) _ 


1 + co 2 



qki(qii m i + q2i m 2) 
|coj 2 - co 2 j + J2£,j cajca 


( 12 . 66 ) 


12.3.7 Two-DOF Response Equations, Frequency Example 



Frequency (Hz) 

Figure 12.20. Two-DOF System Transmissibility 

The system in Figure 12.17 has the parameters previously shown in Table 6.1, with 5% damping for each 
mode. Its natural frequencies are again 20.0 and 42.5 Hz. The transmissibility function for each mass 
relative to the base is given in Figure 12.20 as calculated per equation (12.66). The damping could be 
determined by apply the half-power bandwidth method to the curves if this were measured data. 
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13 Sine Sweep Testing 

13.1 Sine Sweep Introduction 



Time (sec) 

Figure 13.1. Sample Sine Sweep Time History with Very Brief Duration for Educational Purposes Only 

Structures and components may be mounted on shaker tables, similar to those in Figure 4.9 and Figure 
4.10, and subjected to sine sweep base excitation. The essence of a sine sweep test is that the base 
excitation input consists of a single frequency at any given time. The frequency itself, however, is varied 
with time. The sine sweep test may begin at a low frequency and then sweep to a high frequency, or 
vice-versa. Some specifications require several cycles, where one cycle is defined as from low to high 
frequency and then from high back to low frequency. 

The specification might require either a linear or a logarithmic sweep rate. The sweep will spend 
greater time at the lower frequency end if the sweep is logarithmic. The example in Figure 13.1 had a 
logarithmic sweep rate and a constant amplitude. Nevertheless, the specification might require that the 
amplitude vary with frequency. 
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Some potential test objectives are: 

1. Identify natural frequencies and amplification factors or damping ratios 

2. Perform sine sweep before and after random vibration test to determine if any parts 
loosened, etc. 

3. Check for linearity of stiffness and damping by varying the input amplitude 

4. Workmanship screen for defective parts and solder joints 

5. Represent an actual environment such as a rocket motor oscillation 

6. NASA/GSFC typically uses sine sweep vibration for spacecraft testing 


13.2 Sine Sweep Time Domain Equations 


13.2.1 Linear Sweep Rate 

The normalized amplitude function for a linear sweep is 


Y(t) = sin< 7t 




''o 


1 / 


+ f, 


(13.1) 


The test duration is T. The starting and ending frequencies are fj and f 2 , respectively. 

13.2.2 Logarithmic Sweep Rate 

The normalized amplitude function for a logarithmic sweep is 


Y(t) = sin< 


2k 


(fl)(-l + 2 Rt 
Rln(2) 



(13.2) 


The sweep rate in octaves per time is R. Recall the octave formula Section 5.6. A typical sweep rate is 
one octave/minute for sine sweep tests intended for natural frequency identification. 
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13.2.3 Sweep Type Comparison 



Figure 13.2. Sweep Rate Comparison Example 

A comparison of the frequency time variation for two sweep rate types is shown in Figure 13.2. The 
sweeps were performed over the frequency domain from 20 to 2000 Hz over a 400 second duration. 
This is equivalent to a rate of one octave/minute for the logarithmic sweep. 

13.3 Sine Sweep Specification Example 

A vendor has a product that must withstand sinusoidal vibration with an amplitude of 12 G, over the 
frequency domain from 10 Hz to 2000 Hz, with a logarithmic sweep rate. 

The shaker table has a displacement limit of 1.0 inch peak-to-peak, or 0.5 inch zero-to-peak. Recall that 
the displacement limit is a constraint at low frequencies. How should the test be specified? The answer 
is to use a specification with two amplitude segments as shown in Figure 13.3. The first segment is a 
constant displacement ramp, 0.5 inch zero-to-peak. The second segment is a constant acceleration 
plateau at 12 G. 
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Figure 13.3. Sample Sine Sweep Specification 

The cross-over frequency is calculated from the square root of the ratio of the peak acceleration to the 
peak displacement. The respective peaks should be in zero-to-peak format. 


f cross — a / X 


2 n V ^peak/Xpeak 


The cross-over frequency for the specification in Figure 11.2 is 


(13.3) 


f = — 



= 15.3 Hz 


(13.4) 


This is the frequency at which the 12 G acceleration corresponds to a 0.5 in zero-to-peak displacement. 
It is the knee frequency in Figure 13.3. 
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13.4 Sine Sweep Input & Response Example, SDOF System 

The SDOF system in Figure 12.12 has a natural frequency of 40 Hz and amplification factor Q=10. It is 
subjected to a base input of 1 G from 10 to 80 Hz, with a logarithmic rate of 1 octave/min. The input 
and response time histories are shown in Figure 13.4 and Figure 13.15. The peak response is 10 times 
the peak input acceleration. This corresponds to Q=10, but this is the special of a single-degree-of- 
freedom excited into resonance. A better method for damping identification is the half-power 
bandwidth method, which is particularly needed for multi-degree-of-freedom and continuous systems 
due to modal gain factors. 

A transmissibility function can be made relating the response to the input in the frequency domain. This 
is done by tracking the frequency in terms of the zero-crossing rate. The result is shown in Figure 13.6 
along with the half-power bandwidth calculation. Each of the points in Figure 13.6 are approximations 
constrained by the resolution and plot cursor limitations. The points at 37.85 and 41.91 Hz are each 
half-power points, 3 dB less than the resonant peak. The amplification factor is the frequency difference 
between the half-power points divided by the center frequency. 



40.05 Hz 

(41.91-37.85)Hz 


= 9.86*10 


(13.5) 


Recall that Q = 10 corresponds to 5% damping. 
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Figure 13.4. Sine Sweep Input & Response, Total Duration 
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Figure 13.5. Sine Sweep Input & Response, Close-up View 
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Figure 13.6. SDOF Transmissibility Function from Input & Response Time Histories 


13.5 Sine Sweep Input & Response Example, Continuous Beam 



w(t) 

/ / / 7 / 


Figure 13.7. Cantilever Beam, Base Excitation 
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Table 13.1. Cantilever Beam Parameters 


Material 

Aluminum 

Length 

24 inch 

Cross-section 

Circular, 1 inch diameter 

Area Moment of Inertia 

0.04909 in A 4 

Elastic Modulus 

1.0e+07 lbf/in A 2 

Mass/Length 

0.0785 Ibm/in (2.035e-04 Ibf sec A 2/in A 2) 

Amplification Factor 

Q=10 for all modes 


Consider the cantilever beam in Figure 13.7 with the properties in the previous table. 
The fundamental frequency of the cantilever beam using the formula in Table 8.3. 


f, =- 


i 

3.5156 

J(l ,0e + 07 lbf / in 2 ) (o.04909in 4 ) 

2n 

(24 in) 2 

]j 2.035e - 04 lbf sec 2 / in 2 


= 47.7 Hz 


(13.6) 


The first four natural frequencies and their associated parameters are shown in the following table. 


Table 13.2. Cantilever Beam Modal Parameters 


Mode 

Natural 

Frequency(Hz) 

Participation 

Factor 

Modal Mass Ratio 

1 

47.7 

2.99E-03 

0.646 

2 

299 

9.20E-04 

0.198 

3 

837 

3.16E-04 

0.068 

4 

1641 

1.62E-04 

0.035 
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The cantilever beam is subjected to a 1 G sine sweep from 20 to 2000 Hz at a logarithmic rate of one 
octave/minute. The acceleration response is shown in Figure 13.8, as calculated using Reference [20]. 
The corresponding transmissibility function is shown in Figure 13.9. 
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Figure 13.9. Cantilever Beam Transmissibility at Free End 

Again, each of the modes had an amplification factor of 10, but each of the four peaks in Figure 13.9 has 
an amplitude other than 10. The peaks depend on both the amplification factor and the mode shape 
parameters. Nevertheless, the half-power bandwidth method can be used to extract the proper 
amplification values as show for the first and second modes in the following equations. 

The amplification factor for the first bending mode is 

f 47.86Hz 

Af (50.15-45.36)Hz 


The amplification factor for the second bending mode is 


Q = 


Af 


299.9 Hz 
(314.1-284) Hz 


«10 


(13.8) 


The half-power bandwidth calculation could be repeated for the third and fourth modes. 
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13.6 Random Excitation Sources 


Flight Accelerometer Data, Suborbital Launch Vehicle 

n-1-1-1-r 



20 30 40 

Time (sec) 


Figure 13.10. Nonstationary Flight Accelerometer Data 


The liftoff acoustic event occurs during the first two seconds. The vibration then tapers down to a 
lower lever during the subsonic portion of flight. The envelope then increases as the vehicle accelerates 
through the transonic velocity and encounters its maximum dynamic pressure condition in the 10 to 20 
second segment. The level then gradually decreases throughout the remainder of the first stage burn. 
Some transient events occur beginning at 60 seconds due to the attitude control system and other 
sources. 

Consider the following vibration sources: 

1. Turbulent airflow passing over an aircraft wing 

2. Oncoming turbulent wind against a building 

3. Rocket vehicle liftoff acoustics 

4. Earthquake excitation of a building 

One common characteristic of these examples is that the motion varies randomly with time. Thus, the 
amplitude cannot be expressed in terms of a "deterministic" mathematical function. 
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Dave Steinberg wrote [10]: 

The most obvious characteristic of random vibration is that it is nonperiodic. A knowledge of the 
past history of random motion is adequate to predict the probability of occurrence of various 
acceleration and displacement magnitudes, but it is not sufficient to predict the precise magnitude 
at a specific instant. 

Random vibration analysis is thus a matter of probability and statistics. The parameters include: mean, 
standard deviation, histogram, power spectral density, etc. Stationary random vibration is where the 
key statistical parameters remain constant over time. The goal of most shaker table vibration tests is to 
apply a stationary base input to the test unit. This would be the case where the specification was a 
power spectral density. But real-world vibration tends to be nonstationary as show in 
Figure 13.10. 

13.7 Statistical Parameters 

The mean value p of the time history is calculated as 


1 n 


n 


i=l 


where Yj is the instantaneous amplitude and n is the total number of points. 


The variance g z of the time history is 




n 


i=l 


The root-mean-square value Y RMS of the time history is 


(13.9) 


(13.10) 


l rms 


1 n 
1 


b, 1 '' 


i=l 


(13.11) 


A normal distribution is symmetric about its mean value, but some time histories may have skewed 
histograms. The skewness y is calculated as 


Y = 


M_ 

nc 3 


(13.12) 
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Kurtosis is a measure of the "tailedness" of the histogram. The kurtosis k is calculated as 


K 


m-nf 


(13.13) 


Sinusoidal time histories have k = 1.5. Broadband random, stationary, Gaussian time histories have 
K = 3.0. Time histories with higher kurtosis values usually have higher sigma peaks than those 
expected from a Gaussian distribution. 

The crest factor is the ratio of the absolute peak divided by the RMS. 

C — Apeak / Arms ( 13.14) 

The denominator on the right-hand side may be replaced by the standard deviation for the case of zero 
mean. The crest factor for pure sine vibration is always V2. The crest factor for random vibration is 
typically assumed to be 3, but higher crest factors can occur even for closed-loop stationary broadband 
vibration test performed on shaker tables [21]. 

Now consider a single-degree-of-freedom system subjected to a stationary Gaussian broadband random, 
applied force or base excitation. The expected crest factor C e per Reference [22] is 


Note that equation (13.15) can be extended to a multi-degree-of-freedom response by replacing the 
natural frequency with the positive-slope zero-crossing frequency. 

13.8 Random Vibration Histograms 

Again, the histogram divides the time history amplitude points into bins and then counts the points in 
each bin. The fundamental histogram for random vibration is the Gaussian or Normal distribution, also 
known as a bell curve. This function represents the distribution of the instantaneous points for a 
stationary random vibration environment, which may either be broadband or narrowband. White noise 
is a special case of broadband vibration which is present later in this document. The Gaussian histogram 
can be normalized as probability density function such that the total area under the curve is equal to 
one. These functions are useful in vibration analysis for estimating peak responses and fatigue damage. 
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jit - 3o \x-2o n-a |ii |^ + o n + 2a n + 3o 


Figure 13.11. Gaussian Probability Density Function. 
(https://www.kdnuggets.com/ 2018/ 08/ basic-statistics-python-probability.html) 

Again, the mean value is [J and the standard deviation is Q. The mean value in random vibration is 
usually zero, either because the vibrating object has zero net displacement or because the 
accelerometer signal has been AC coupled. 

The Gaussian distribution equation is 


P(x) = 



f 

exp 

V 


z m!1 

2a 2 


(13.16) 


The probability that the variable x will be between any two given limits is obtained by integrating the 
area under the curve, which requires the gamma function. Probability examples are shown in the 
following tables. 


Table 13.3. Gaussian Distribution, Probability of Occurrence Between Limits 


Statement 

Probability 

Ratio 

Percent 

-a < x < +a 

0.6827 

68.27% 

-2a < x < +2a 

0.9545 

95.45% 

-3a < x < +3a 

0.9973 

99.73% 
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Table 13.4. Gaussian Distribution, Probability of Occurrence Outside Limits 


Statement 

Probability 

Ratio 

Percent 

1x1X7 

0.3173 

31.73% 

| x | > 2a 

0.0455 

4.55% 

1 X | > 3(7 

0.0027 

0.27% 



Figure 13.12. Rayleigh Distribution 


The Rayleigh distribution is a continuous probability distribution for positive-valued random variable 
with the equation 


p(A) = — 7 exp 

CT 


^-A 2 ^ 


A > 0 


(13.17) 
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The X-axis represents a scale factor times the standard deviation, as denoted by A which is used to 
denote the absolute value of x. 

The area under the curve from zero to a given A value represents probability of occurrence. The total 
area under the curve is one. The probability of exceeding the A value is one minus the area as shown 
in 13.5. The probability that a given peak will be within zero and la is 39.35%. The probability of 
exceeding la is 60.65%. 

An SDOF response time history can be reduced to a subset of points representing local peaks. A peak 
occurs at a point where the initial slope has the opposite polarity of the ending slope. A peak can be 
negative or positive regardless of whether it is a local minimum or maximum. The peak distribution will 
tend to be a Rayleigh distribution for random excitation, whether by applied force or base input. 
Characterizing the distribution is important for fatigue analysis. 


Table 13.5. Rayleigh Distribution, Probability of Occurrence Outside Limits 


X 

Prob [ A > Xo ] 

0.5 

88.25 % 

1.0 

60.65 % 

1.5 

32.47 % 

2.0 

13.53 % 

2.5 

4.39% 

3.0 

l.l 1 % 

3.5 

0.22 % 

4.0 

0.034 % 


Thus, 1.11 % of the peaks will be above 30 for a signal whose peaks follow the Rayleigh distribution. 
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13.9 White & Pink Noise 



Figure 13.13. White Light Passing through Prism 


White light passed through a prism produces a spectrum of colors. White noise is a type of analogous 
random vibration. White noise and pink noise are two special cases of random vibration. 

White noise is a random signal which has a constant power spectrum for a constant frequency 
bandwidth. It is thus analogous to white light, which is composed of a continuous spectrum of colors. 
Static noise over a non-operating TV or radio station channel tends to be white noise. In addition, 
Commercial white noise generators are designed to produce soothing random noise which masks 
household noise as a sleep aid. 

Pink noise is a random signal which has a constant power spectrum for each octave band. This noise is 
called pink because the low frequency or “red” end of the spectrum is emphasized. Pink noise is used in 
acoustics to measure the frequency response of an audio system in a particular room. It can thus be 
used to calibrate an analog graphic equalizer. Waterfalls and oceans waves may generate pink-like noise. 
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13.10 White Noise Input & Response Example 



o 
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o 
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0 10 20 30 40 50 60 


Acceleration Response fn=600 Hz Q=10 



Time(sec) 


Response Histogram 

100 

50 

o 

0 

-50 

-100 
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Figure 13.14. White Noise Base Input Time History and SDOF Response with Histograms 


Stationary white noise base input excitation is applied to the SDOF model in Figure 12.12. The system 
has a natural frequency of 600 Hz and an amplification factor of Q equal to 10, or 5% damping. The 
input time history is shown in Figure 13.14, with a Gaussian distribution and a standard deviation of 10 
G. The response is calculated using the Smallwood ramp invariant digital recursive filtering relationship 
[17]. The response likewise has a Gaussian distribution. 
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Figure 13.15. SDOF Response to White Noise Base Input, Close-up View 
Table 13.6. SDOF Response to White Noise Statistics 


Parameter 

Value 

Mean 

0G 

Std Dev 

21.8 G 

RMS 

21.8 G 

Skewness 

0 


Parameter 

Value 

Kurtosis 

3.0 

Crest 

Factor 

4.28 

Maximum 

93.3 G 

Minimum 

-90.6 G 


A brief segment from the input and response are given in Figure 13.15. The white noise input is 
broadband random. The response tends to vibrate at 600 Hz natural frequency, but the amplitude 
envelope varies randomly with time. The response is narrowband random. The response statistics for 
the entire duration are given in Table 13.6. 

The response mean is zero. The RMS and standard deviation values are thus equal per equation (17) 
and have a value of 21.8 G. The highest input absolute peak is 93.3 G. The ratio of the peak relative to 
the standard deviation is called the “crest factor,” which is 4.28 in this example. This value can also be 
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represented as 4.28CJ. A common assumption is that the crest factor is 3G, which is three times the 
standard deviation. But 4Q and 5Q peaks are common in stationary random vibration. 


Local Peak Absolute Values Peak Histogram 
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60 
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Figure 13.16. SDOF Response to White Noise, Absolute Peak Values with Histogram 


The histogram of the absolute local peaks has a Rayleigh distribution. The expected crest factor for the 
sample problem is 


C= ^2 In [(600 Hz) (60 sec)] 


_ 0.5772 _ 

^2 In [(600 Hz)(60 sec)] 


4.71 


(13.18) 


The crest factor for sample problem was 4.28, which is 10% lower than the expected value. 

13.11 Crest Factor Study 

As a numerical experiment, the white noise synthesis and response calculation steps for the previous 
example were repeated 20,000 times. A crest factor was calculated for each trial. A histogram of the 
crest factor is shown in Figure 13.17. The histogram has a shape somewhat similar to a Rayleigh 
distribution, but further investigation is needed to properly characterize it. The crest factors ranged 
from 4.14 to 6.62, with a mean of 4.76. This experimental mean value was 1.1% higher than the 
theoretical 4.71 value. Again, stationary random vibration time history with Gaussian distributions 
typically have peaks ranging from 4a to 5a, but higher crest values may also occur. This is true both for 
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broadband inputs and narrowband responses. Designing structures and components for 3CT vibration 
loads is thus insufficient, even though this is a common rule-of-thumb. 



4 4.5 5 5.5 6 6.5 7 

Crest Factor 

Figure 13.17. Flistogram of Response Peaks for 20,000 Trials 


13.12 Non-Gaussian Random Vibration 

Stationary Gaussian random time histories can be generated in closed-loop shaker table tests. But time 
histories measured in the field usually depart from the Gaussian ideal. A theoretical example is shown in 
Figure 13.18. Nevertheless, the response of an SDOF system to a non-Gaussian random input will have 
some tendency to be Gaussian due to the “central limit theorem.” 

The central limit theorem in probability theory establishes that, in some situations, when independent 
random variables are added, their properly normalized sum tends toward a normal distribution even if 
the original variables themselves are not normally distributed. 
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Figure 13.18. Non-Gaussian Broadband Random Time History & Response 


The base input has standard deviation=1.0, kurtosis=9.3 and crest factor=ll.l. The response statistics 
are shown in the following table. 


Table 13.7. SDOF Response to Non-Gaussian Time History 


Parameter 

Value 

Kurtosis 

3.5 

Crest 

Factor 

6.9 

Maximum 

12.7 G 

Minimum 

-14.6 G 


Parameter 

Value 

Mean 

0 G 

Std Dev 

2.12 G 

RMS 

2.I2G 

Skewness 

0 
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13.13 HALT/HASS Testing with Non-Gaussian Random Vibration 



Figure 13.19. HALT/ HASS Testing of Circuit Boards, Combined Thermal Cycling and Vibration 
(https://www.dfrsolutions.com/ hubfs/ Resources/ Sherlock/ How-Sherlock-Can-Help-You-Now-and- 

Future.pdf) 


Highly Accelerated Life Testing (HALT) is a process for ruggedization of preproduction products. 
Highly Accelerated Stress Screening (HASS) is the production screen for products once they have been 
characterized in HALT, typically performed on 100% of units per lot. The random vibration 
environment is nonstationary, non-Gaussian vibration with high kurtosis. The overall GRMS of the 
vibration is controlled, but the spectral content is uncontrolled. The thermal cycling and vibration may 
be applied separately or simultaneously depending on the test phase. Note that HALT/HASS is the 
offspring of the classical Environmental Stress Screening (ESS) process. 
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Figure 13.20. HALT/ HASS Chamber Vibration Impactors 
(http://www.elitetest.com/ sites/ default/files/ downloads/ mechanical_shock_and_vibration_testing.pdf) 

The underside of the chamber is connected to hammer impactors which create a repetitive shock time 
history. 



Time (sec) 


Counts x 10 5 


Figure 13.21. HALT Time History Simulation with Histogram 
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A time history was synthesized to represent a typical HALT chamber vibration environment, as shown 
in Figure 13.21. 

The kurtosis is 8.9 and the crest factor is 9.3. The time history appears stationary, but the close-up view 
in Figure 13.22 reveals that it is composed of a series of transient pulses. 


HALT Simulation 5 GRMS 



Figure 13.22. HALT Time History Simulation, Close-up View 
The vibration is a series of repetitive shock pulses due to the hammer strikes underneath the table. 
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Figure 13.23. Electronics Life Cycle Bathtub Curve (Wikimedia) 

A successfully screening process for mass-produced electronics will precipitate latent defects to 
eliminate the “infant mortality” failures shown in Figure 13.23. The defects could be bad solder joints, 
loose connectors or screws, bent clamps, flawed pieces parts, etc. An example of a bad solder joint is 
shown in Figure 13.24. 



The purpose of subjecting electronic 
components to HALT/HASS is to 
uncover flaws such as bad solder 
joints. The electronic components 
can be powered and monitored 
before, during, and after the tests to 
verify that they are functional. 


Figure 13.24. Cracked Solder J oints in Piece Part with "J Leads" 
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14 Sample Rate & Aliasing 

14.1 Data Collection Accuracy 

The wagon-wheel effect, or stroboscopic effect, is an optical illusion in which a spoked wheel appears to 
rotate at a different speed than its true rotation. The wheel can appear to rotate more slowly than the 
true rotation, appear stationary, or even appear to rotate in the opposite direction from the true 
rotation. The wagon-wheel effect is seen in film or television depictions of stagecoaches or wagons 
in Western movies, but the effect may also be seen in videos of helicopter rotors and aircraft propellers. 

This optical illusion is a form of aliasing. It can occur in a video recording where the frame rate is too 
slow relative to the rotational or oscillation frequency of the moving object. A similar problem can 
occur when accelerometer data is digitized unless precautions are taken. 

Engineers collect accelerometer data in a variety of settings both in the field and in lab tests. The 
accelerometers measure the data in analog form. The accelerometer may have an integral mechanical 
lowpass filter. Furthermore, the signal conditioning unit may have an analog lowpass filter. Lowpass 
filtering of the analog signal is necessary to prevent aliasing. Eventually, the data is passed through an 
analog-to-digital converter. 

The proper lowpass frequency and sampling rate must be selected to ensure that the digitized data is 
accurate. There are two rules-of-thumb which must each be met. 

14.2 First Requirement 

The first requirement is that the sampling rate must be greater than the maximum analysis frequency. 
Industry guidelines are summarized in the following table, where 

(minimum sampling rate) > ( N )( maximum analysis frequency ) (14.1) 


Table 14.1. Sample Rate First Requirement 


Analysis Type 

N 

Frequency Domain 

2 

Time Domain 

10 


The frequency domain requirement is based on the need for at least two time-domain coordinates per 
cycle to resolve a sine wave. The frequency domain analysis thus extends up to the Nyquist frequency 
which is one-half the sample rate. Note that some conservative references specify an N of 2.5 for 
frequency domain calculation. 
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A sampling rate of 100 kHz is thus required for a shock response spectrum (SRS) analysis extending to 
10 kHz. The SRS function is covered in Section 0. Recall that the SRS is calculated in the time domain. 

The IES Handbook for Dynamic Data Acquisition and Analysis [23] gives the following guidelines: 

Unlike other spectral quantities evolving from the discrete Fourier transform computations, 
the SRS is essentially a time domain quantity. Hence, the digital sampling rate given by Rs=l/At 
introduces errors beyond those associated with aliasing about the Nyquist frequency. Thus, 

Rs must be high enough to accurately describe the response of the SRS oscillators. To 
minimize potential error, it is recommended that the SRS computations be performed with a 
sampling rate of Rs > 10 fh, where fh is the highest natural frequency of the SRS computation. 

14.3 Second Requirement 

Shannon’s sampling theorem states that a sampled time signal must not contain components at 
frequencies above the Nyquist frequency. 

Again, the Nyquist frequency is equal to one-half the sampling rate. Shannon's theorem applies to 
frequency domain analysis. 

Thus, the second requirement is that the sampling rate must be greater than the maximum frequency 
present in the source energy at the measurement location. This requirement is independent of the 
maximum analysis frequency. 

The guidelines are summarized in the following table, where 

(minimum sampling rate) > ( M )( maximum frequency in source energy ) (14.2) 

Table 14.2. Sample Rate Second Requirement 


Analysis Type 

M 

Frequency Domain 

2 

Time Domain 

10 


Note the similarity between Table 14.1 and Table 14.2. 
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14.4 Anti-aliasing Filter 



Figure 14.1. Data Acquisition Flow Chart 

The challenge is that the maximum frequency in the source energy maybe unknown, which is the usual 
case for field data. The solution is to use an analog lowpass filter prior to digitization. This filter is also 
referred to as an anti-aliasing filter. The filter cut-off frequency is typically set at, or slightly above, the 
maximum analysis frequency. 

The IES Handbook for Dynamic Data Acquisition and Analysis gives the following guidelines. 

Let 


f c be the cutoff frequency 
fN be the Nyquist frequency 

• A lowpass anti-aliasing filter with a cutoff rate of at least 60 dB/octave should be used for the 
analog-to-digital conversion of all dynamic data. 

• With a 60 dB/octave cutoff rate, the half-power point cutoff frequency of the filter should be set 
at f c < 0.6 fn. 

• If the anti-aliasing filter has a more rapid cutoff rate, a higher cutoff frequency can be used, but 
the bound f c < 0.8 fn should never be exceeded. 
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Figure 14.2. Pyrotechnic Shock Test of an Avionics Component 

Mild detonating cord is source shock and is coiled on the underside of the plate. Accelerometers are 
mounted on the plate to record the shock level. The component is powered and monitored during the 
test. The source shock may have energy components exceeding 100 kHz. Use of an analog anti-aliasing 
filter is critical. Otherwise, high frequency energy will be erroneously folded down to lower frequencies 
in the shock response spectrum. 

14.5 Aliasing Example 


Sine Function 20 Hz 



Figure 14.3. A 20 Hz Sine Function Sampled at 200 Hz 
There are 10 points per period. The Nyquist frequency is 100 Hz. 
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Sine Function 180 Hz 



Figure 14.4. An 1800 Hz Sine Function Samples at 2000 Hz 

The red markers show that there are 1.11 points per 1800 Hz cycle. The Nyquist frequency is 1000 Hz. 
Aliasing occurs! 



Figure 14.5. Two Signals Digitized at 2000 Hz 
The two signals both appear as 200 Hz although only one truly is. 
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Fourier Magnitude 180 Hz Signal Aliasing Error 
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Figure 14.6. Fourier Magnitude of Each Signal Sampled at 200 Hz 


The 180 Hz signal is folded about the Nyquist frequency which is 100 Hz. The resulting energy is 
deposited at 20 Hz. Fourier transforms are covered in Section 15. 

As an example, consider a pyrotechnic shock test with a time domain maximum analysis frequency of 10 
kHz. The sample rate should be at least 100 kHz. The Nyquist frequency is 50 Hz. Assume a lowpass 
filter with a 60 dB/octave cutoff rate. The maximum cutoff frequency should thus be 30 kHz. These 
are values to satisfy basic requirements. A better approach would be to set the sample rate at 250 kHz 
or more if the measurement is going to be in the “near field” relative to the source. This higher sample 
rate can be readily achieved in a ground test, but not in flight where the telemetry bandwidth is limited. 

As another example, consider a vibration test with a maximum analysis frequency of 2 kHz in the time 
domain. The sample rate should be at least 20 kHz. The Nyquist frequency is 10 kHz. The maximum 
cutoff frequency should thus be 6 kHz. 
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15 Fourier Transforms 

15.1 Fourier Transform Equations 

15.2 Textbook Fourier Transform 

The Fourier transform is a method for identifying the spectral content of a time history. This transform 
has already been used previously in this document, including for the channel beam vibration response in 
Figure 11.3 and Figure 11.4. This section will give background equations starting with the rarefied 
mathematician approach and finishing with the practical engineering version. 

The textbook Fourier transform X(f) for a continuous time series x(t) is defined as 


X(f) = J x(t)exp(-j27if t)dt 


(15.1) 


Note that j = ^ -1 . And both the time and frequency domains are considered to extend from minus 
to plus infinity in their respective spheres. 

There is also an inverse Fourier transform defined as 

x(t) = J_ oo X(f)exp(+j27ift)df (15.2) 

As an aside, the corresponding Euler’s equation is 

exp(± j2nf t) = cos(27if t)± jsin( 2nf t) (15.3) 

There are several impractical quirks begat from this transform and inverse pair. One is that a signal as 
simple as a pure sine function will have a spectral component at a negative of its true frequency, as 
shown in Figure 15.1. How is an engineer to design and test for a negative frequency? In addition, 
measured time histories have a distinct start and stop time with a finite duration. Next, integration is 
analogous to multiplication. Assume that x(t) represents an acceleration time history. The resulting 
Fourier transform from equation (15.1) would return a transform with amplitude dimension 
[acceleration-time]. What would this be? Velocity? 

The last peculiarity is that measured time histories are strictly real. The resulting Fourier transform of 
a measured signal from equation (15.1) would be complex, with real and imaginary components. Now 
take this Fourier transform and insert it into the inverse formula in (15.2). The resulting time history 
will be complex. Fortunately, this concern is easily resolved because the imaginary components will all 
be zero. 
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Imaginary X(f) 



f=0 

Figure 15.1. Mathematicians' Fourier Transform of a Sine Function 
The time history behind the imagery Fourier transform equation (15.1) is 


x(t) = Asin^27tf 


(15.4) 


The real Fourier transform is omitted because it is zero for a pure sine function with zero phase angle. 
The symbol 8 represents the Dirac delta function, a mathematical oddity unsuitable for engineering 
design and test purposes. 

Note the following traits for any given signal: 

1. The real Fourier transform is symmetric about the f=0 line. 

2. The imaginary Fourier transform is anti-symmetric about the f=0 line. 


15.2.1 Discrete Fourier Transform 

The discrete Fourier transform rectifies some, but not all, impracticalities of the textbook Fourier 
transform. Engineers collect data using data acquisition systems which return digitized outputs for 
measured analog inputs. The integral in equation (15.1) thus needs to be exchanged for a series, which 
is easily done. The resulting discrete Fourier transform is 


N-l 


V 1 J f • 2 71 , 

Xu = — > < x n exp -j—nk 

k N n t 0 l \ J N 


, for k = 0,1,...,N-l 


(15.5) 
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Its corresponding inverse is 


N-l 

Z 

n=0 


‘n = Z ”1 X k ex pf+J^ nk 


for n = 0,1,..., N-l 


(15.6) 


Eureka! Equation (15.5) does away with negative frequencies. The starting frequency will now be zero. 
Also, the Fourier transform amplitude dimension will now be the same that of the time history. The 
time domain dimension of acceleration will thus carry directly over to the Fourier transform. And both 
the time and frequency domain functions may now be defined over finite domains. Alas, this discrete 
pair will need further tweaking before the equations are ready for practical use, but progress has been 
made. 

Recall the beloved sine function from equation (15.4). Set its amplitude at 1 G and its frequency at 1 Hz. 
Digitize it with a sample rate of 32 samples per second with a 512 second duration. The frequency step 
in the Fourier transform will be (1/512) Hz. The imaginary discrete Fourier transform is shown in the 
following figure. 


Imaginary DFT y(t) = 1 sin [ 2 n (1 Hz) t ] 
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Figure 15.2. Discrete Fourier Transform of a Sine Function 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


139 










enDAQ 


Chapter: Fourier Transforms 


The sine function has a frequency of 1 Hz. Again, the Fourier transform of a pure sine function is 
imaginary and antisymmetric. The antisymmetric reflection occurs about 16 Hz which is the Nyquist 
frequency, one-half the sample rate. Does the discrete Fourier transform imply that that a hypothetical 
component subjected to this environment needs to be designed and tested for 31 Hz as well as 1 Hz? 
How would the peak at 31 Hz change if the sample rate were increased to 64 samples per second? Is 
this some manifestation of Heisenberg uncertainty principle? And why is the absolute amplitude at 1 Hz 
equal to 0.5 but the sine function amplitude is 1? 

The discrete Fourier transform in equation (15.5) will not be abandoned. It reappears in certain 
calculations such as going back-and-forth between impulse response functions and transfer functions, 
convolution applications, inverse Fourier transforms, etc. These intermediate concepts will be covered 
in a future document. But a more practical transform is needed presently. 

15.2.2 Useful Discrete Fourier Transform 

Sidestepping the quandaries raised in the previous section, let us modify the discrete Fourier transform 
into a practical tool. Henceforth, the Fourier transform will only be shown up to the Nyquist frequency. 
And the Fourier transform will be represented in terms of magnitude rather than real and imaginary 
components. The phase may also be shown with magnitude if needed. And the amplitudes at all 
nonzero frequencies will be doubled. A more elegant of way of expressing this is that the energy 
components above the Nyquist frequency will be folded about the Nyquist frequency to lower 
frequencies. The energy at 31 Hz in 15.2 would thus be deposited at 1 Hz. Note that the spectral 
magnitude at zero frequency is a special case since it remains the simple mean value of the signal. 

The full amplitude Fourier transform up to the Nyquist frequency is 




, N—1 

— Z x n> for k = 0 

iN n=0 


2 _ 

N 


V 1 f - 2n v 

2^i x n exp -j— nk 

n=() ^ 


for k = 1, ...,—-1 
2 


(15.7) 
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Frequency (Hz) 

Figure 15.3. Full Amplitude, Discrete Fourier Transform of a Sine Function, Up to Nyquist Frequency 

The beloved sine function’s Fourier magnitude is plotted in the above figure using equation (15.6). A 
clean spectral peak of 1 G appears at 1 Hz. The phase angle is omitted because it would only be valid at 
the sinusoid’s own frequency. It would be noise at other frequencies. There is a need for phase angle in 
other situations, such as modal testing, comparing the lag from input to response, etc. 

15.3 Leakage & Hanning Window 

15.3.1 Leakage Error 

There are numerous error sources associated with the Fourier transform, including "leakage,” which is a 
smearing of energy throughout the frequency domain. Leakage results when both of the following 
conditions are present: 

1. The signal is taken over a finite duration 

2. The signal is "non-periodic" in the time record 

Both these conditions are usually present in engineering data, so leakage is common. For example, 
leakage occurs if a Fourier transform is calculated for a non-integral number of sine function cycles. 
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Figure 15.4. Sine Function and Memory Buffer Example 1 


Consider that a data acquisition system is used to monitor a sine function over a 10 second duration. 
The sine function has an amplitude of 1 G and a frequency of 1 Hz, as shown in the top subplot. The 
sample rate is 32 samples per second. Now assume that the data acquisition system has a limited 
memory buffer and measures only three cycles as shown in the middle subplot. Note that the time 
history amplitude is zero at the start and end of the record, with an integer number of cycles. In 
essence, the Fourier transform will correctly assume that the original signal is a series of three-cycle 
segments as shown in the bottom subplot, such that no leakage occurs as shown in Figure 15.6. 
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Figure 15.5. Sine Function and Memory Buffer Example 2 


Repeat the previous example but with a tighter memory buffer that only measures 2 1/2 cycles as shown 
in the middle subplot. Leakage will occur due to the non-integer number of cycles. The Fourier 
transform will assume that the original signal is a distorted series of these segments as shown in the 
bottom subplot, as shown in Figure 15.6. 
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Figure 15.6. Fourier Magnitude Comparison, Both Examples 


The top subplot shows a clean spectral line at 1 Hz for the case of three whole cycles. The bottom 
subplot shows that the energy is smeared across numerous bands as a result of leakage error from the 
non-integer number of cycles. 

On a related note, the top subplot has a spectral line available at exactly 1 Hz. But the bottom subplot 
lacks a line at 1 Hz and instead has lines at 0.8 and 1.2 Hz. Recall that the frequency step is the inverse 
of the duration. The three-cycle spectral increment is thus 0.333 Hz, whereas the two and half cycle 
step is 0.4 Hz. 
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15.3.2 Hanning Window 


Hanning Window w(t) 



Time 

Figure 15.7. Hanning Window 


The leakage error can be reduced by subjecting the time history to an appropriate window. Two 
common types of windows are the rectangular window and the Hanning window. 

The rectangular, or flat, window leaves the time history data unmodified and is equivalent to no window 
at all. A rectangular window is appropriate for transient data or nonstationary data, including the type 
of shock time histories presented in Section 0. It does not reduce leakage. 

One of the most common windows is the Hanning window, or the cosine squared window. It works 
best with stationary data. This window tapers the time history data so that the amplitude envelope 
decreases to zero at both the beginning and end of the time segment. The Hanning window w(t) can be 
defined the is plotted in Figure 15.7. 


w(t) = 1 — cos^ 




,0<t<T 


(15.8) 


A normalization factor of the square root of (8/3) may be applied to the Hanned data to compensate for 
lost energy, with the goal of maintaining the original overall RMS level. 
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Figure 15.8. Sample Sine Function, with Non-Integer Number of Cycles 


A sine function before and after Hanning window application is show in the above figure. The non¬ 
integer number of cycles will cause a leakage error which is mitigated by the Hanning window although 
not eliminated. The Fourier transforms of each are shown in Figure 15.9. 
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Fourier Magnitude, y(t) = 1.0 sin ( 2n t), Af = 0.105 Hz 



Figure 15.9. Fourier Spectra, Sine Function Before and After Hanning Window 


Ideally, the Fourier transform would have a single, discrete line at 1 Hz with an amplitude of 1 G. A 
challenge is that there are spectral lines available at 0.945 Hz and 1.05 Hz but none at 1.0 Hz. Both the 
rectangular and Hanning Fourier transforms have some leakage error as a result. The rectangular 
window produces more leakage error than the Hanning window. Thus, the Hanning window is 
recommended for stationary data. 

15.4 Fast Fourier Transform 


15.4.1 FFT Characteristics 

Table 15.1. Suitable Record Lengths for the FFT Calculation 
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The discrete Fourier transform requires a tremendous amount of calculations for digital signals with 
large numbers of points. A time history with M coordinates would require M 2 complex multiplication 
steps. The discrete Fourier transform can be carried out by a Fast Fourier transform (FFT) method, 
however. The method is based on a time series with a number of points equal to 2 N , where N is an 

integer. The FFT requires Mlog2M complex multiplication steps, where M = 2 N . Some length examples 
are shown in Table 15.1. 


15.4.2 FFT Example 

Now consider a time history with 1,000,000 points. A regular Fourier transform would require 10 12 
complex multiplication steps. On the other hand, an FFT would only require approximately 2(10 7 ) 
steps. Thus, the FFT achieves the calculation in 1/50,000 th of the time. 

The above example is not quite correct. Again, the FFT is based on a time series with 2 N coordinates. 
Note that 2 19 = 524,288 and 2 20 = 1,048,576. But a time history with 1,000,000 points falls between 
these two cases. 

There are some options for dealing with a time history that is not an integer power of 2. One option is 
to truncate the time history. This should be acceptable if the data is stationary. In the above example, 
the time history would thus be truncated to 524,288 points. The second option is to pad the time 
history with trailing zeroes to bring its length to an integer power of 2. A problem with this option is 
that it artificially reduces the amplitude of the Fourier transform spectral lines. A third option is to 

carefully perform linear interpolation on the data so that its length becomes 2 N . Or the problem can 

be circumvented altogether by configuring the data acquisition system to collect exactly 2 N points. This 
was a feature of some vintage systems. 


15.4.3 FFT Algorithm Background 

The Fourier transform requires multiplying the time points by complex weighting functions. These 
functions have a certain symmetry for the case of time series with 2 N points. The calculation can thus 
be simplified by taking advantage of this symmetry using concepts such as the unit circle and the 

butterfly algorithm. A few introductory steps are given here. Define a weighting factor W n . 


W n = exp 


-J- 


. 2nn 


N 


(15.9) 


The discrete Fourier transform becomes 


X(k) 


isjxoow*}. for k = 0,l 


N-l 


(15.10) 
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The matrix representation is 


X(0) 


w° 

w° 

w° • 

.. W 0(N-1) 

x(0) 


X(1) 

_ 1 
” N 

w° 

w 1 

w 2 

.. W KN-1) 

x(l) 


X(2) 

w° 

w 2 

w 4 • 

.. W 2(N-1) 

x(2) 

(15.11) 

X(N-1)_ 


w° 

W KN-1) 

W 2(N-1) . 

.. W (N-1)(N-1) 

_x(N-l)_ 



Note that the weighting matrix in equation (15.11) is symmetric. Also note 

W° = 1 


The Fourier transform equation in matrix form simplifies to 


X(0) 


X(1) 

1 

X(2) 

_X(N —1)_ 

N 


1 

1 

W 1 

w 2 

w 2 

w 4 

W KN-1) 

W 2(N-1) 


1 

W KN-1) 
W 2(N—1) 


x(0) 

x(l) 

x (2) 


W (N-1)(N-1) 


x(N-l) 


The weights can be arranged geometrically as show in the following figure. 


(15.12) 


(15.13) 



Figure 15.10. Unit Circle for N=8, Polar Coordinates with Real and Imaginary Axes 
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The units circle shows some emerging symmetries. The real part of W 1 is equal to that of W 7 . The 

imaginary part of W 1 is equal to that of W 3 . Additional symmetries can be identified by noting that 
some values are the negative of others. Further development of the FFT algorithm is beyond the scope 
of this document but is given in References [24] and [25]. 

15.5 Electrical Transformer Hum 



Figure 15.11. Transformer Utility Box 

Magnetic fields from alternating current induce vibration in certain materials, such as iron, steel, and 
ferrite alloys. This effect is called magnetostriction. Transformer hum is a common example of 
magnetostriction. This hum is a potential nuisance, and it is an energy loss mechanism. 
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Figure 15.12. Step-up Transformer 
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A basic transformer has two sets of windings, the primary and the secondary. Certain transformers may 
have secondary windings or taps. The primary connects to the power source. The secondary connects 
to the load usually at a lower voltage. The windings are wrapped around a core, which may be iron or 
ferrite, depending on the design requirements. An iron core may be laminated with layers of iron and 
non-conducting material stacked together. A laminated iron core is used to reduce the creations of 
eddy currents in the iron core that would dissipate the energy being transferred from the primary coil 
to the secondary coil in the form of heat. 

An alternating current is input to the primary winding. The current produces a magnetic flux, or 
magnetic lines of force in the core. The resulting magnetic field induces an alternating current in the 
secondary winding. The number of turns on each winding determines the output voltage from the 
transformer. The output voltage from the secondary is proportional to the ratio of the turns on the 
windings. 



Figure 15.13. Magnetostriction Cycle 


The magnetic field stresses the core material. The stress causes the core volume to expand regardless 
of the polarity. The core contracts to its equilibrium volume as the magnetic flux returns to zero. The 
material expands and contracts twice per each magnetic cycle, as shown in Figure 15.25. The core thus 
vibrates at 120 Hz in response to the changing magnetic field from 60 Hz AC power. The deformation 
amplitude may be so small that it is imperceptible to the unaided human eye. Nevertheless, it may be 
large enough to generate a hum. The hum frequency is 120 Hz, with integer harmonics. A hypothetical 
strain gage mounted on the core would output a rectified version of the curve in Figure 15.25. 
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Figure 15.14. Transformer Flum Fourier Spectra 

The Fourier magnitude from a sound output of the transformer in Figure 15.11 is shown in Figure 15.14, 
with its 120 Hz peak and integer harmonics. 

15.6 Apache Helicopter Flyover Data 
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Figure 15.15. Apache AH-64 Helicopter 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


152 
























enDAQ 


Chapter: Fourier Transforms 


The Apache is a twin-engine army attack helicopter, as shown in Figure 15.15. It entered service with 
the US Army in 1984. 

Helicopter noise consists of a complex mixture of sounds. The repeating impulse noise from the rotor 
blades is the dominant source of distinct spectral tones. The fundamental tone occurs at a characteristic 
frequency depending on the number of blades and their rotation speed. The blade-passing frequency 
(BPF) is the shaft rotation frequency times the number of rotor blades. Noise occurs at the BPF as well 
as at integer harmonics. Both the main rotor and the tail rotor generate this type of noise. 

Rotor-vortex interaction noise is another source. This source is also referred to as blade-vortex 
interaction (BVI) noise. This source is significant in descent or level flight at low and medium velocities. 
Vortices form in the wake of the blades. A given rotor blade can run into the tip vortex shed by a 
preceding blade. This causes a distinctive type of noise: the annoying "blade-slap" of helicopters with 
slow-turning rotors, or a sharp fluttering noise for rotorcraft with fast-turning rotors. 

Furthermore, the velocity at the tips of the advancing blades in high-speed forward flight may approach 
Mach 1. Shock waves begin to form at this condition, leading to shock induced flow separation. The 
shock waves and turbulent flow generate further noise. This source is referred to as highspeed impulsive 
(HSI) noise. 

In addition, most helicopters are powered by turbine engines. Engine and gearbox noise consist of 
distinct spectral peaks as well as broadband random noise. 



Figure 15.16. Apache Flyover Acoustic Time History 

A segment of an acoustic time history from an Apache flyover is shown in the above figure. The signal is 
a summation of sine tones, mostly from the main and tail rotor blade passing frequencies. 
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Figure 15.17. Apache Flyover Acoustic Fourier Spectra with Main Rotor BPFs Identified 

The measured blade passing frequency is 21 Hz with integer harmonics. The Apache’s main rotor has 
four blades. The apparent main hub frequency is thus 5.25 Hz. The actual main hub frequency is = 4.84 
Hz per published data in Reference [26]. The apparent frequency is higher due to a Doppler shift with 
the Apache flying toward the microphone. 


Fourier Magnitude Apache Flyover 



Figure 15.18. Apache Flyover Acoustic Fourier Spectra with Tail Rotor BPFs Identified 
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The Apache tail rotor has four blades. The blades, however, are not oriented 90° (perpendicular) from 
each other as in most helicopters. Specifically, one set in front of the other at a 55° angle. The 
supplementary angle is 125°. This unusual arrangement is required because the two sets of blades use a 
"Delta-Hinge" which allows the blades to simultaneously flap and feather. The four blades appear to 
behave as two for the tail rotor blade pass frequency. 

The purpose of this design is to render the blades more efficient and to reduce noise. The main rotor 
hub and the tail rotor hub are driven by the same source. Thus, there is a fixed ratio between the spin 
rates of each rotor. The tail rotor spin rate is about 4.86 times higher than the main rotor rate per 
Reference [26]. Again, the measured, Doppler-shifted, main rotor spin rate is 5.25 Hz, or 315 rpm. The 
tail rotor spin rate is 25.5 Hz, or 1530 rpm, per the fixed ratio. The unique configuration of the blades 
yields a 2X blade passing frequency at 51 Hz. The blades behave as two pairs rather than as four 
individual blades. Tail rotor integer harmonics of this frequency occur at 102, 153, and 204 Hz. Recall 
that 153 Hz was the dominant peak in the Fourier transform in Figure 15.17, due to the tail rotor’s 
interaction with the main rotor wake. 

15.7 Turboprop Engine 



Figure 15.19. Bombardier Q400 Aircraft 


The author once flew from Spokane to Seattle, Washington on an Alaska/Horizon Air Bombardier Q400 
aircraft, similar to the one in Figure 15.19. He made an audio recording using his Android phone. The 
Fourier magnitude is shown in Figure 15.20. This model aircraft has two Pratt & Whitney Canada 
PW150A turboprop engines. The PW150A engine/propeller rotation rate during takeoff and climb is 
1020 RPM but is throttled back at cruise altitude to 850 RPM, or 14.17 Hz. There are six blades on each 
engine, so the blade passing frequency is 85 Hz. This frequency and its integer harmonics are clearly 
visible in the spectral analysis. 
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Figure 15.20. Turboprop Engine Noise Fourier Spectra 


15.8 Waterfall FFT 

15.8.1 Introduction 

A waterfall FFT is a 3D plot with axes of time, frequency, and amplitude. It is useful for visualizing how 
spectral peak amplitudes and frequencies change with time. A spectrogram is a related 2D plot that 
shows this same data but with amplitude represented in terms of contour colors. 
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15.8.2 Solid Rocket Motor Example 



BsaafiEF 


Figure 15.21. Solid Rocket Motor 

Solid rocket motors may have standing pressure oscillations which form in the combustion chamber 
cavities, caused by vortex-shedding and other effects. This condition is called “Resonant Burn” or 
“Thrust Oscillation.” The sinusoidal oscillation frequency may sweep downward as the cavity volume 
increases due to the conversion of propellant to exhaust gas. A typical upper-stage motor susceptible 
to this problem in shown in Figure 15.21. Flight data from an accelerometer mounted on a bulkhead 
adjacent to this motor during its resonance burn is shown in Figure 15.22. The sinusoidal character is 
apparent in the close-up view. 
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Solid Rocket Motor, Resonant Burn, Flight Data 
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Figure 15.22. Solid Motor Resonance Burn, Accelerometer Data, Full and Close-up Views 



Solid Rocket Motor, Resonant Burn, Flight Data 
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Waterfall FFT 



Figure 15.23. Solid Motor Resonance Burn, Accelerometer Data, Waterfall FFT 

The resonant burn spectral peaks begin at 520 Hz and then sweep downward to 450 Hz. The 
corresponding amplitude varies with both frequency and time probably due to structural resonance 
effects as measured by the bulkhead accelerometer. The causes of peaks at 280 and 350 Hz are 
unknown. 


15.8.3 Wind Chime Acoustics 
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Figure 15.24. Wind Chime, Unsealed Sound Pressure, Waterfall FFT 
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Recall the wind chime example from Section 8.1.7. The bending frequencies for the tallest chime were 
238, 657, 1288 and 2127 Hz. A waterfall of the chime’s sounding is shown in Figure 15.24. The 
fundamental tone is very persistent with low damping. The second mode makes a very small 
contribution. The third mode initially has the highest generated sound, but the tone decays quickly. 
The fourth mode’s spectral peaks are barely visible in the plot. 

15.9 White Noise Fourier Transform 

The Fourier transform is an excellent tool for resolving both the frequency and amplitude of a pure 
sinusoidal signal, aside from leakage concerns. Alas, it is a poor choice for random vibration. This is the 
motivation for representing random vibration in terms of power spectral density as presented in Section 
16. The problem with the Fourier transform approach for random vibration relates to the dependence 
of the amplitude on the frequency resolution as shown in the following example. Consider a stationary 
white noise signal. An engineer may choose the entire signal or a representative segment for the 
Fourier transform calculation. A longer segment duration gives a finer frequency resolution. There 
would thus be more points in the Fourier transform magnitude to represent the signal’s energy. The 
corresponding amplitude of each point would be, on average, lower than if a shorter segment were used 
with a wider frequency resolution and fewer points. 
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Figure 15.25. White Noise Time History 

A stationary white noise signal is shown in Figure 15.25. The standard deviation of the entire time 
history is 1 G. The standard deviation of the segment from 0 to 1 second is also 1 G, approximately. 
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Separate Fourier transforms are taken from the entire duration and for the 0 to 1 second segment, as 
shown in Figure 15.26. 



0 10 20 30 40 50 60 70 80 90 100 

Frequency (Hz) 

Figure 15.26. White Noise Fourier Transform 


Each curve represents a standard deviation of 1 G, equivalent to 1 GRMS given zero mean. The overall 
GRMS value can be calculated by taking the “square root of the sum of the squares” from each curve. 
The blue curve has fewer points, by a factor of one-tenth, to carry the energy than the red curve. The 
blue curve’s average amplitude is thus higher than the denser red curve. 
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16 Power Spectral Density 

16.1 Decibels & Slopes 

16.1.1 Differences 

Power spectral density functions are sometimes specified in terms of decibels. The dB difference 
between two levels A & B in units of G 2 or G 2 /Hz is 

AdB = 10 log (A / B) M6 h 


A 6 dB increase would quadruple a G 2 /Hz level. A 6 dB decrease would reduce the G 2 /Hz level by a 
ratio of one-fourth. 

The dB difference between two levels C & D in units of G or GRMS is 


AdB = 20 log (C / D) ( , 62) 

A 6 dB increase would double a G level. A 6 dB decrease would reduce the G level by a ratio of one- 
half. 


16.1.2 Slopes 

Consider a frequency domain plot with a straight-line segment in log-log format with coordinates (fi, yi) 
and (f 2 , y2). 

The slope N between the coordinate pair is 


N = x °g[y2/yi] 
log[ f 2 / f i] 


(16.3) 


The line equation for this pair is 
y<f) = y,(f N /f, N ) 

The dB/octave slope is for units of G 2 or G 2 /Hz is 


AdB / oct = 10 N log 2 


(16.4) 


(16.5) 
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The dB/octave slope is for units of G or GRMS is 

AdB / oct = 20 N log 2 


(16.6) 


16.2 PSD Specifications 



Figure 16.1. Navmat P-9492 Acceleration PSD 


The power spectral density in Figure 16.1 is taken from Reference [27]. The lowercase g convention is 
carried over from this source, but uppercase G is used elsewhere in the present document This level is 
intended as a base input for a shaker table test to screen components for parts and workmanship 
defects. It is representative of similar levels found in various military and NASA standards. The Y-axis 
unit is more properly grms 2 /Hz. The specification consists of three straight line segments in log-log 
format. The overall grms value is the “square root of the area under the curve.” The integration 
approach is to take the grms 2 area under each of the three segments. Then sum the areas and take the 
square root. This requires special integration formulas due to the log-log format which depend on the 
slope from equation (16.3). 
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The area ai for a segment is 


a i = < 
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(r \ 
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v f i. 


, for N = -1 


( 16 . 7 ) 


The overall level for m total segments is then 


a RMS“^X a i ( 16 . 8 ) 

Performing these steps for the Navmat P-9492 specification in Figure 16.1 requires some preliminary 
work, because the specification did not explicitly give the first and last amplitudes. The first coordinate 
amplitude is 0.01 g 2 /Hz by inspection of the graph. Identifying last coordinate amplitude requires a two- 
step calculation. The slope is 


N = 


-3dB 

101og2 


- 1.0 


The final coordinate amplitude is 


y = 


0.04 g 2 / Hz 
(350 Hz) -1 


(2000 Hz) 1 = 0.007 G 2 /Hz 


( 16 . 9 ) 


( 16 . 10 ) 


The final amplitude is whimsically referred to as “Bond.” The overall grms calculation is summarized in 
the following table. 
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Table 16.1. Navmat P-9492 Overall Level Calculation 


Frequency 

Domain 

Slope N 

grms 2 

20 to 80 

i 

1.5 

80 to 350 

0 

10.8 

350 to 2000 

-1 

24.4 


Total 

36.7 


The overall level is the square root of 36.7 grms 2 which is 6.06 grms. 

16.3 PSD Integration in the Frequency Domain 

An acceleration PSD has a corresponding velocity PSD and a displacement PSD. Let 

APSD = Acceleration PSD 
VPSD = Velocity PSD 
DPSD = Displacement PSD 

The integration formulas are 

VPSD((d) = APSD(oo) / co 2 (16.11) 

DPSD(co) = APSD(co) / ca 4 (16.12) 

The resulting velocity and displacement PSDs for the Navmat specification are shown in Figure 16.2 and 
Figure 16.3, respectively. These plots are useful to determine whether a given acceleration PSD test can 
be performed on a certain shaker table given the table’s velocity and displacement limits. 
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Frequency (Hz) 


Figure 16.2. Navmat P-9492 Velocity PSD 



Frequency (Hz) 

Figure 16.3. Navmat P-9492 Displacement PSD 
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16.4 PSD Calculation Methods 

Power spectral density functions of measured data may be calculated via three methods: 

1. Measuring the RMS value of the amplitude in successive frequency bands, where the signal 
in each band has been bandpass filtered 

2. Taking the Fourier transform of the autocorrelation function 

3. Taking the limit of the Fourier transform X(f) times its complex conjugate divided by its 
period T as the period approaches infinity. 

The filtering option is most instructive method and is covered in Section 20. The autocorrelation option 
is the Wierner-Khintchine approach and is omitted from this document. A modified form of the third 
option is covered in this section and is the preferred method for computational efficiency. 

The one-sided power spectral density function Xpsp)(f) is calculated from the discrete Fourier 
transform X(f) as 


Xps D (0 = 


lim 

Af —» 0 


1 X(f)X*(f) 

2 Af 


(16.13) 


The Fourier transform in extends from zero frequency to the Nyquist frequency, which is one-half the 
sample rate. Each Fourier transform has peak dimension. The one-half factor is needed to convert the 
amplitude from peak 2 /Hz to rms 2 /Hz. 

The frequency step is finite in practice and is the inverse of the total measured duration. 

Af =1/T (16.14) 

This frequency step is the smallest sine wave frequency which can be resolved, where this sine wave has 
one period equal to the total duration. The frequency step is linked to the number of statistical- 
degrees-of-freedom. The reliability of the power spectral density data is proportional to the number of 
degrees. A wider Af gives greater PSD confidence in terms of smoothing the spectral components. 

The number of degrees N^f is defined as 


N dof =2BT (16.15) 

The variable B is the bandwidth of an ideal rectangular filter. It is equal to the frequency step Af for an 
ideal rectangular filter. The product BT is unity and the number of degrees-of-freedom is always equal 
to two for a given record assuming a rectangular filter. This number may be increased for a given 
record by subdividing it into a set of small segments, where each segment has two degrees. The total 
number of degrees is then equal to twice the number of segments. This subdivision approach tends to 
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smooth the PSD. The penalty, however, is that the frequency resolution widens as the record is 
subdivided which could smear narrow peaks. These tradeoffs are shown by example in Table 16.2. 


Table 16.2. Time History with 4096 Samples Taken over 16 Seconds, Rectangular Filter 


Number of 
Records 

NR 

Number of 
Time Samples 
per Record 

Period of 
Each Record 
Ti 

(sec) 

Frequency 

Resolution 

Bj = 1 / Tj 
(Hz) 

dof per 
Record 

= 2Bj Ti 

Total dof 

i 

4096 

16 

0.0625 

2 

2 

2 

2048 

8 

0.125 

2 

4 

4 

1024 

4 

0.25 

2 

8 

8 

512 

2 

0.5 

2 

16 

16 

256 

1 

1 

2 

32 

32 

128 

0.5 

2 

2 

64 

64 

64 

0.25 

4 

2 

128 


Additional rows could be added to the above table but are omitted for brevity. Selecting the optimum 
processing strategy for a given time history is a matter of engineering judgment. Continuing with the 
table example, a frequency step of 4 Hz with 128 dof might be a good choice for demonstrating that the 
broadband base input for a shaker table test complied with the specified PSD’s tolerance bands. On the 
other hand, a narrow frequency resolution from the first or second row may be needed if the purpose 
of analyzing the PSD it to identify lightly-damped, low-frequency resonant peaks. 

Recall that a rectangular filter leaves the time history data unmodified. Leakage error may be a concern, 
however, as discussed in Section 15.3. A Hanning window can be applied to each data segment to 
reduce leakage, but it is a non-ideal filter which reduces the number of degrees-of-freedom. The key to 
applying the Hanning window is to subdivide the data into overlapping segments. Nearly 90% of the 
degrees-of-freedom are recovered with a 50% overlap, as shown in Figure 16.4. 
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Original Sequence 



Segments, 
Hanning Window, 
50% Overlap 


Figure 16.4. White Noise and Subdivision into Overlapped Segments with Hanning Windows 


16.5 PSD Calculation Examples 
16.5.1 White Noise PSD 

Recall the white noise time history in Figure 15.25. A PSD is calculated for this signal with Hanning 
window and 50% overlap for each of two frequency resolution cases. The PSDs are shown in Figure 
16.5. 
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Frequency (Hz) 
Figure 16.5. White Noise PSDs 


The blue curve has a finer frequency resolution and more variation than the red curve. The red curve is 
smoother due to its wider frequency resolution and higher dof number. The red curve effectively 
averages the PSD points into wider bands. Note that pure white noise would have a flat PSD curve. 
The time history in Figure 13.5 departs from this ideal because it was taken over a short, 10 second 
duration. It was also band-limited via low-pass filtering such that it has a roll-off beginning at about 85 
Hz. Filtering is covered in Section 20. 
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16.5.2 Taurus Auto PSD 



Figure 16.6. 2003 Ford Taurus Auto 

The author mounted an accelerometer in his auto’s console as shown in Figure 16.6. The vehicle was 
driven on a highway at 65 mph. An accelerometer time history from this road test is given in Figure 
16.7. The PSD is shown in Figure 16.8, with a narrow frequency resolution to aid in spectral peak 
identification. The fundamental mode damping is calculated via the half-power bandwidth method in 
Figure 16.9. 



Figure 16.7. Taurus Time History 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


171 













































Accel (G z /Hz) 


eilDAQ 


Chapter: Power Spectral Density 



Figure 16.8. Taurus PSD, Af =0.3 Hz, 6 dof 



Figure 16.9. Taurus PSD, Half-Power Bandwidth 
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Q = 


_f_ 

Af 


1.5 Hz 

(1.9-0.88) Hz 


«1.5 


(16.16) 


The Taurus auto spring-mass Frequency is 1.5 Hz with an amplification factor equivalent to 33% 
damping, which is typical damping for shock absorbers. Common automobile natural frequencies are 
given in the following table. 


Table 16.3. Automobile Spring-Mass Frequencies 


Vehicle 

Fundamental Frequency 

Passenger Car 

1 to 1.5 Hz 

Sports Car 

2 to 2.5 Hz 

Hummer 

4.5 Hz 


The spectral peaks in Figure 16.8 at 14.6 and 29.1 Hz can be explained as tire imbalance frequencies per 
the following calculations steps. 

a. Assume 25 inch tire outer diameter at 65 mph 

b. Circumference = TT (25 inch) = 78.5 inch 

c. 65 mph =1144 in/sec 

d. Speed/Circumference = (1144 in/sec) / 78.5 in = 14.6 Hz 

e. 2X harmonic = 29.1 Hz 

16.6 PSD Time History Synthesis 


Table 16.4. Time History Synthesis Steps to Meet PSD Specification 


Step 

Description 

i 

Generate a white noise time history 

2 

Take the FFT 

3 

Scale the FFT amplitude per the PSD for each frequency 

4 

The time history is the inverse FFT 

5 

Use integration, polynomial trend removal, fade in and out, and differentiation so 
that corresponding mean velocity and mean displacement are both zero 

6 

Scale the time history so that its GRMS value matches the specification’s overall 

GRMS value 

7 

Take a PSD of the synthesized time history to verify that it matches the PSD 
specification 
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The method for synthesizing a time history to satisfy a PSD are shown in Table 16.4 . These steps are 
used to synthesize an acceleration time history to satisfy the Navmat P-9492 specification. The resulting 
broadband random time history is shown in Figure 16.10, with its normal distribution histogram. It 
began as white noise but was modified such that the final time history is no longer white noise. Its final, 
corresponding PSD is shaped and defined over a finite frequency domain, as shown in Figure 16.11. 
White noise would have a flat PSD in contrast. 


Acceleration Time History Synthesis 

40 

30 

20 

10 

o 

~o 0 
o 

4 

-10 

-20 

-30 

-40 

0 20 40 60 80 100 120 140 160 180 0 5 10 



Histogram 



T\me(sec) Counts x 1G 5 

Figure 16.10. Acceleration Time History and Histogram 
Table 16.5. Acceleration Time History Statistics 


Parameter 

Value 

Mean 

0 G 

Std Dev 

6.06 G 

RMS 

6.06 G 

Skewness 

0.011 


Parameter 

Value 

Kurtosis 

3.01 

Crest 

Factor 

5.375 

Maximum 

32.56 G 

Minimum 

-30.93 G 
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Figure 16.11. PSD Verification 

The synthesized time history satisfies the specification well within the tolerance bands. The frequency 
step is 2.84 Hz with 1024 statistical degrees-of-freedom. 
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Figure 16.12. Synthesized Time History, Acceleration, Velocity and Displacement 

Each of the three response time histories has a stable oscillation about its zero baseline. Each also has a 
brief fade in and out, which could be seen in a close-up view of the start and finish. 

Now assume an SDOF system with (fn=200 Hz, Q=10) subjected to the Navmat P-9492 PSD base input, 
using the model in Figure 12.12. Solve for the acceleration response in the time domain using the 
synthesized base input. The numerical engine is the Smallwood ramp invariant digital recursive 
relationship. The response time history is shown in Figure 16.13. 
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Figure 16.13. SDOF Response to Synthesized Base Input 
Table 16.6. Acceleration Response Time History Statistics 


Parameter 

Value 

Kurtosis 

3.02 

Crest 

Factor 

4.548 

Maximum 

51.36 G 

Minimum 

-50.95 G 


Parameter 

Value 

Mean 

0 G 

Std Dev 

11.3 G 

RMS 

11.3 G 

Skewness 

0.003 


The overall response of 11.3 GRMS is very close to the theoretical values from the frequency domain 
methods which will be shown in Sections 16.7 and 16.8. The theoretical crest factor for this case per 
equation (13.18) is 4.71, slightly higher than the 4.548 value in the table. 
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The time domain synthesis method could whimsically be referred to as a Rube Goldberg approach, after 
the famous inventor and cartoonist. 


Base Input 




Figure 16.14. Input and Response Acceleration including Response Velocity, Close-up View 

The base input is broadband random vibration. The response acceleration and velocity time histories 
are narrowband random. The SDOF system prefers to oscillate at its natural frequency. The positive 
slope zero cross rate is 199.6 Hz. 
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Figure 16.15. Input and Response PSDsfor Navmat P-9492 Synthesis 

The response PSD tracks the input at the low frequency end with nearly unity gain. The resonant 
response occurs at and near the 200 Hz natural frequency. The energy above V2 times the natural 
frequency is attenuated. Compare Figure 16.15 with the same set of PSD curves derived from the 
frequency domain method which will be shown in Figure 16.18. 
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Figure 16.16. Power Transmissibility for Input and Response PSDs 

The peak power transmissibility reaches nearly 100 G 2 / G 2 at 200 Hz which is equivalent to Q 2 where 
Q=10. But this is a special case of an SDOF system subjected to base excitation. A more robust 
method for estimating the Q value, as if this were experimental data, is to use the half-power bandwidth 
method shown in Figure 16.17. 
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Frequency (Hz) 

Figure 16.17. Power Transmissibility with Half-Power Bandwidth Points 
The half-power bandwidth method yields the following amplification factor. 


Q = 


f_ 

Af 


199.4 Hz 

(209.6-188.9)Hz 


= 9.63 * 10 


( 16 . 17 ) 


16.7 Miles Equation for Base Excitation 

Consider an SDOF system subjected to base excitation where the input PSD is white noise over the 
frequency domain from 0 to infinity Hz, using the model in Figure 12.12. The resulting acceleration 
X G rms can determined by Miles equation from Reference [10]. The constant power spectral 
density amplitude is represented by A with unit of (G A 2/Hz). 


*GRMS ( f n> Q) 


^ 71 ^ 




f n QA 


An equivalent form is 


XGRMs( f n> S) = . 


71 


v^y 




f n A 


( 16 . 18 ) 


( 16 . 19 ) 
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The Miles equation is widely used due to its simplicity, but its assumption of white noise over an infinite 
domain does not exist in physical reality. A rule-of-thumb states that it may be used if the input PSD is 


flat within one-octave on either side of the natural frequency. In practice, this rule is used with some 


compromise. 

Assume an SDOF system with (fn=200 Hz, Q=10) subjected to the Navmat P-9492 PSD base input. The 
overall response is 



= 11.2 GRMS 


(16.20) 


A more robust method for performing the response calculation is given in Section 16.8. 

16.8 General Method for Base Excitation 

Real-world PSD specifications are shaped and have lower and upper frequency limits. Miles’ equation 
cannot account for a PSD with ramps and plateaus. The problem is exacerbated if the PSD is from 
narrowband measured data with peaks and dips. These practical characteristics can be readily 
accounted for by applying the power transmissibility function to the base input PSD and then by doing a 
point-by-point multiplication calculate the response PSD. The overall response is then the square root 
of the area under the response PSD curve. 

Recall the SDOF transmissibility function from equation (12.47). 



(16.21) 


The power transmissibility is equal the transmissibility squared. The response PSD is calculated from 
the power transmissibility to the input PSD. 



(16.22) 


*GRMS ( f n - Q) = 



(l-Pi 2 ) 2 +(2^ Pi ) 2 


l + (^Pj) 2 


Y PSD( f i) Af i - Pi= f i/ f n (16.23) 
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Equation (16.23) is cumbersome but is readily implement via a software program. It does not appear to 
have a name in the literature but is referred to as the “general method” in this document. 



Frequency (Hz) 

Figure 16.18. General Method SDOF Response to Navmat P-9492 Base Input PSD 


Repeat the example of an SDOF system with (fn=200 Hz, Q=10) subjected to the Navmat P-9492 PSD 
base input using the general method. The overall response in Figure 16.18 agrees with the Miles result 
in equation. It is also very close to the time domain overall response in Table 16.6. The general method 
would yield a more accurate result if the natural frequency fell on either of the input PSD slopes. 

Compare Figure 16.18 with the same set of PSD curves derived from the time domain synthesis in 
Figure 16.15. 
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17 Vibration Response Spectrum 

17.1 VRS Basics 

The vibration response spectrum overall level is calculated using equation (16.23) by varying the natural 
frequency for a given Q value. It can be plotted in terms of the overall 1-sigma, 3-sigma, and peak values 
as a function of natural frequency. The GRMS and 1-sigma values are taken to be equal because the 
mean is assumed to be zero. The expected peak value, or crest factor, is calculated from equation 
(13.15) and is the only of three response metrics that depends on duration. 

The vibration response spectrum is demonstrated in the following example. First calculate the SDOF 
response to the Navmat P-9492 PSD base input for natural frequencies at 100, 200, and 300 Hz, with 
Q=10, using the model in Figure 12.12. 



Figure 17.1. SDOF Response to Navmat P-9492 Base Input PSD for Three Natural Frequencies 
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Natural Frequency (Hz) 


Figure 17.2. Vibration Response Spectrum of Navmat P-9492 Base Input PSD 


The calculation of the response PSDs in Figure 17.1 could be repeated for a larger set of natural 
frequencies. This process can be automated via a software program. The resulting VRS is shown in 
Figure 17.2. 

The Figure 17.1 legend’s three coordinates are shown explicitly in the plot. A duration of 180 seconds 
is assumed, but this only affects the peak curve. The peak curve shows that peaks up to 4 and 5-sigma 
may occur over the duration. 
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17.2 VRS Comparison of Two PSDs 



Figure 17.3. PSD Comparison 

The following example is actual case history. A launch vehicle avionics component was tested to the red 
curve in the above figure. The blue curve represents the flight data from an accelerometer mounted 
adjacent to the component. The flight data has three peaks that exceeded the test level. Was the 
component under-tested? The natural frequency of the component is unknown but could have been at 
any of the three exceedance peaks. Assume that each curve represents the same duration for simplicity. 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


186 


























enDAQ 


Chapter: Vibration Response Spectrum 



Figure 17.4. Vibration Response Spectrum Comparison 


The test level was 8.6 dB higher that the flight data in terms of overall level. The vibration response 
spectrum comparison shows that the component would have a higher response to the test level than to 
the flight level across the entire natural frequency spectrum from 20 to 2000 Hz. The component test 
level enveloped the flight data with margin to spare. 

17.3 PSD Envelope Derivation via VRS 

The vibration response spectrum can also be used an indirect method for deriving a PSD envelope for a 
measured environment. The PSD envelope will have a few coordinates in log-log format as a 
simplification. The goal is to derive the least PSD possible that will still envelop the measured PSD in 
terms of the respective VRS curves. This allows for some peak clipping as demonstrated in the 
following example for an actual case history. There are some industry rules-of-thumb for peak clipping, 
but the VRS is a more thorough approach. 
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Frequency (Hz) 

Figure 17.5. Suborbital Launch Vehicle Flight PSD and Envelope 


The flight data from a launch vehicle internal bulkhead location is shown in the above figure. The 
spectral peaks were unexpected. A realization was later made that a small access cover had been 
prematurely jettisoned at liftoff. The result was that cylindrical module housing the bulkhead became a 
large whistle. This was an anomaly which could be corrected for the next flight in this vehicle series. 
Nevertheless, it presented an opportunity to develop the enveloping method. The following example 
proceeds assuming the whistle condition would be repeated. The proposed envelope is shown with the 
flight data in Figure 17.5. The VRS justification for the envelope is given in Figure 17.6. Assume the 
flight and envelope duration are the same for simplicity. The envelope derivation process uses trial-and- 
error with random number generation, along with some convergence routines. 
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Figure 17.6. Vibration Response Spectrum Comparison 
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18 Shock Testing & Analysis 

18.1 Shock Isolator Photo Gallery 



Figure 18.1. Titan II Missile Silo, Launch Control Room 

The control room is mounted underground via huge isolation springs. A typical spring is shown in the 
background. The purpose is to isolate the control room from mechanical shock and vibration in the 
event of a nuclear strike above the launch site. The springs allow 18 inches of relative displacement. The 
control room could thus carry out a retaliatory strike, as ordered by the U.S. president. This site is 
located south of Tucson, Arizona. It has been decommissioned and is now a museum. 
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18.2 Simple Drop Shock 



Figure 18.2. iPhone 6 Accidental Drop 

The first iPhone 6 was purchased in Perth, Australia on September 19, 2014. The event was covered by 
a live TV report. The buyer mishandled the phone as he unboxed it. The phone survived the drop 
onto the ground but may have had some fatigue or fracture damage. 

Portable electronic devices (PEDs) are expected to survive multiple drops. Most original equipment 
suppliers specifying between 30 and 50 drops. Recommended test methods are given in Reference [28]. 



/////// 



x 


Figure 18.3. Drop Shock Analytical Model 

The accidental drop shock of a component onto a hard surface is difficult to model accurately. The item 
may undergo rigid-body translation and rotation in each of its six degrees-of-freedom during freefall. 
The item may have a nonlinear response with plastic deformation during impact. It may or may not 
bounce. Furthermore, a box-shaped object may strike the ground on any of its corners, edges or faces. 

A very simple method, as a first approximation, is to assume that the object is a linear, undamped, 
single-degree-of-freedom subjected to initial velocity excitation as it strikes the ground and remains 
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attached to it via its spring. The object then undergoes free vibration in his configuration. The initial 
velocity is calculated using a familiar physics formula where the change in kinetic energy is equal to the 
change in potential energy due to gravity. 

Assume that the object is dropped from rest. The initial velocity as it strikes the ground is 
x(0) = 72 gAh 
The equation of motion is 
mx + kx = 0 
Or equivalently 

2 

X + 03 n X = 0 

The peak displacement is 

| x(t )| = 7H“ 

co n 

The peak velocity is equal to the initial velocity in equation (18.1). 

The peak acceleration is 

|*(t)| = o) n V2gAh (18.5) 


(18.1) 


(18.2) 


(18.3) 


(18.4) 


An example is shown in the following table for three natural frequency cases. 


Table 18.1. Peak Response Values for 36 inch Drop Height 


Natural Freq 
(Hz) 

Displacement (in) 

Velocity (in/sec) 

Acceleration (G) 

200 

0.133 

167 

543 

600 

0.044 

167 

1630 

1000 

0.027 

167 

2710 
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18.3 Classical Shock 


Half Sine 



Terminal Sawtooth 



Time (sec) 


Versed Sine 




Time (sec) 


Figure 18.4. Idealized Classical Pulse Examples 

Classical pulses are the simplest base excitation pulses. They are deterministic and can be represented 
by simple mathematical functions. They are typically one-sided. An SDOF system’s response to a 
classical pulse can be solved for exactly using Laplace transforms. 

Four classical pulse types are shown in Figure 18.4. Other types include rectangular and trapezoidal 
pulses. These pulses do not necessarily represent real field environments, but they are still used 
throughout industry to test equipment ruggedness for convenience. 

Shock tests are performed on military equipment [26] to: 

a. Provide a degree of confidence that materiel can physically and functionally withstand the 
relatively infrequent, non-repetitive shocks encountered in handling, transportation, and service 
environments. This may include an assessment of the overall materiel system integrity for safety 
purposes in any one or of the handling, transportation, and service environments. 

b. Determine the materiel's fragility level, in order that packaging may be designed to protect the 
materiel's physical and functional integrity. 

c. Test the strength of devices that attach materiel to platforms that can crash. 
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Potential equipment failure modes due to shock excitation include: 

a. Materiel failure resulting from increased or decreased friction between parts, or general 
interference between parts. 

b. Changes in materiel dielectric strength, loss of insulation resistance, variations in magnetic and 
electrostatic field strength. 

c. Materiel electronic circuit card malfunction, electronic circuit card damage, and electronic 
connector failure. (On occasion, circuit card contaminants having the potential to cause short 
circuit may be dislodged under materiel response to shock.) 

d. Permanent mechanical deformation of the materiel resulting from overstress of materiel 
structural and nonstructural members. 

e. Collapse of mechanical elements of the materiel resulting from the ultimate strength of the 
component being exceeded. 

f. Accelerated fatiguing of materials (low cycle fatigue). 

g. Potential piezoelectric activity of materials, and materiel failure resulting from cracks in 
fracturing crystals, ceramics, epoxies, or glass envelopes. 



Figure 18.5. Drop Shock Test Machine, Initial Velocity Excitation 
(Courtesy of Lansmont) 
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Classical pulse shock testing has traditionally been performed on a drop tower. The component is 
mounted on a platform which is raised to a certain height. The platform is then released and travels 
downward to the base, which has pneumatic pistons to control the impact of the platform against the 
base. In addition, the platform and base both have cushions for the model shown. The pulse type, 
amplitude, and duration are determined by the initial height, cushions, and the pressure in the pistons. 
This is a textbook example of case where the initial potential energy of the raised platform and test item 
are converted to kinetic energy. The final velocity of the freefall becomes the initial velocity of the 
shock excitation. 


Acceleration 



Velocity 



Figure 18.6. A 50 G, 11 msec, Terminal Sawtooth Pulse for Shaker Shock Test 

Classical pulse shock testing can sometimes be performed on shaker tables, but there some constraints. 
The net velocity and net displacement must each be zero. Also, the acceleration, velocity and 
displacement peaks must each be within the shaker table stroke limits. Pre and post pulses are added to 
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classical pulses to meet these requirements. A hypothetical terminal sawtooth pulse suitable for shaker 
shock testing is shown in Figure 18.6. 

18.4 Half-Sine Shock Example 

Consider a single-of-freedom is subjected to a 50 G, 11 msec half-sine pulse applied as base excitation 
per Figure 12.12. Set the amplification factor to Q=10. Allow the natural frequency to be an 
independent variable. Solve for the absolute response acceleration. The equation of motion is for the 
relative displacement z is 


.. _ 9 fAsin(7Tt/T) , fort<T 

z + 2^<o n z + <o n z = | ' 0 fo ; t>T (18.6) 

The primary response occurs during the half-sine pulse input. The residual response occurs during the 
quiet period thereafter. The total response is the combination of primary and residual. The exact 
response for a given time can be calculated via a Laplace transform. Note that the quiet period solution 
is free vibration with initial velocity and displacement excitation. 



Figure 18.7. Fourier Transform Half-Sine Shock 

A common misunderstanding is to regard the half-sine shock pulse as having a discrete frequency which 
would be the case if it were extended to a full-sine pulse. The Fourier transform in Figure 18.7 shows 
that the half-sine pulse has a continuum of spectral content beginning at zero and then rolling-off as the 
frequency increases. There are also certain frequencies where the magnitude drops to zero. The 
magnitude represents the acceleration, but the absolute magnitude depends on the total duration 
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include the quiet period after the pulse is finished. The post-pulse duration was 10 seconds in the above 
example, for a total of 10.011 seconds. 


50 
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Acceleration (fn=10 Hz, Q=10) 
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Figure 18.8. SDOF System Response to Half-Sine Shock, 10 Hz 

The response in the above figure has an absolute peak value less than the peak input. This is an isolation 
case. The response positive and negative peaks occur after the base input pulse is over. 

Acceleration (fn=75 Hz, Q=10) 



Figure 18.9. SDOF System Response to Half-Sine Shock, 75 Hz 
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The response in Figure 18.9 has an absolute peak value that is 1.65 times the peak input. This is 
resonant amplification case. The absolute peak response occurs during the base input pulse. 


Acceleration (fn=400 Hz, Q=10) 



Figure 18.10. SDOF System Response to Half-Sine Shock, 400 Hz 

The response converges to the base input as the natural frequency becomes increasingly high. This 
becomes a unity gain case. The system is considered as hard-mounted. 

The peak results from the three cases are shown in Table 18.2. The calculation can be repeated for a 
family of natural frequencies. The peak acceleration results can then be plotted as a shock response 
spectrum (SRS) as shown in Figure 18.11. The peak relative displacement values can likewise be plotted 
as an SRS as shown in Figure 18.12. 
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Figure 18.11. Acceleration SRS, 50 G, 11 msec Half-Sine Pulse 
The two curves in Figure 18.11 contain the coordinates in Table 18.2. 

The initial slope of each SRS curve is 6 dB/octave indicating a constant velocity line. The curves also 
indicate that the peak response can be lowered by decreasing the natural frequency. A low natural 
frequency could be achieved for a piece of equipment by mounting it via soft, elastomeric isolator 
bushings or grommets. But a lower natural frequency leads to a higher relative displacement as shown 
in Figure 18.12. 


Table 18.2. Summary of Peak Response, 50 G, 11 msec, Half-Sine Base Input 


Natural 

Frequency (Hz) 

Peak Positive (G) 

Absolute Value of 
Peak Negative (G) 

10 

20.3 

17.3 

75 

82.5 

65 

400 

51.8 

4.38 
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Figure 18.12. Relative Displacement SRS, 50 G, 11 msec Half-Sine Pulse 


The curves in the above figure could be used for designing isolator mounts for a component. The 
isolators must be able to take up the relative displacement without bottoming or topping out. There 
must also be enough clearance and sway space around the component. 
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Figure 18.13. Acceleration SRS, 50 G, 11 msec Terminal Sawtooth Pulse 


The positive and negative SRS curves are reasonably close for the terminal sawtooth pulse shown in 
Figure 18.13. In contrast, the positive and negative SRS curves for the half-sine pulse in Figure 18.11 
diverge as the natural frequency increases above 80 Hz. The terminal sawtooth pulse is thus usually 
preferred over the half-sine pulse for classical shock testing. 

Another means of visualizing the SRS concept is given in Figure 18.14. 
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Figure 18.14. Half-Sine Shock Applied as Base Input to Independent SDOF Systems 

The systems are arranged in order of ascending natural frequency from left to right and subjected to a 
common half-sine base input. The Soft-mounted system on the left has high spring relative deflection, 
but its mass remains nearly stationary. The Hard-mounted system on the right has low spring relative 
deflection, and its mass tracks the input with near unity gain. The Middle system ultimately has high 
deflection for both its mass and spring. The peak positive and negative responses of each system are 
plotted as a function of natural frequency in the shock response spectrum. 
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18.5 Response to arbitrary excitation 


Recall from Section 12.3 that the response of a single-degree-of-freedom system to base excitation can 
be expressed in terms of a second order, ordinary differential equation for the relative displacement Z 

for a base acceleration y. 



(18.7) 


Equation (18.7) can be solved via Laplace transforms if the base acceleration is deterministic such as a 


half-sine pulse. A convolution integral is needed if the excitation varies arbitrary with time. The 
convolution integral is computationally inefficient, however. An alternative is to use the Smallwood 
ramp invariant digital recursive filtering relationship [17], [19]. 

Again, the recursive filtering algorithm is fast and is the numerical engine used in almost all shock 
response spectrum software. It is also accurate assuming that the data has a sufficiently high sample rate 
and is free from aliasing. One limitation is that it requires a constant time step. 

The equation for the absolute acceleration is 



ex p[ - 2£,co n At]x 




(18.8) 


+ 


exp(-2£,co n At)-— 

v j 



The damped natural frequency is 



(18.9) 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


203 














enDAQ 


Chapter: Shock Testing & Analysis 


The digital recursive filtering relationship for relative displacement is omitted for brevity but is available 
in Reference [19]. 

The relationship in equation (18.8) is recursive because the response at the current time depends on 
the two previous responses, which are the first two terms on the righthand side of the equation. This is 
a feedback loop in terms of control theory. The relationship is filtering because the energy at and near 
the natural frequency is amplified whereas higher frequency energy above V2 times the natural frequency 
is attenuated. See Figure 12.14. 

18.5.1 Shock Response Spectrum Test Specification Objective 


NASA Mid-Field Shock Data 



Time (milliseconds) 

Figure 18.15. Mid-Field Pyrotechnic Shock Time History 

Consider the measured mid-field shock time history in the above figure as taken from Reference [29]. 
The time history would be essentially impossible to reproduce in a test lab given that it is a high- 
frequency, high-amplitude complex, oscillation pulse. The aerospace practice instead is to derive an SRS 
to represent the damage potential of the shock event. The test conductor may then use an alternate 
pulse to satisfy the SRS specification within reasonable tolerance bands. This is an indirect method of 
achieving the test goal. There are some limitations to this approach. One is that the test item is 
assumed to be linear. Another is that it behaves as a single-degree-of-freedom system. Nevertheless, 
this method is used in aerospace, military and earthquake engineering fields, for both analysis and testing. 
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Natural Frequency (Hz) 


Figure 18.16. Mid-Field Shock Response Spectrum & Envelope 


The customary approach is to draw a conservative envelope over the measured SRS. The ramp-plateau 
format is the most common, although there are variations. The enveloping process shown in the above 
vibration is very conservative in the mid frequency domain. An additional dB uncertainty factor may be 
needed to develop the envelope into a test specification, given that the SRS envelope is derived from a 
single time history. Industry standards, such as Reference [30], give guidelines for the dB factor. 
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18.6 Seismic Shock 

18.6.1 Seismic Waveforms 



Figure 18.17. Mother Earth 

The Earth experiences seismic vibration. The fundamental natural frequency of the Earth is 309.286 
micro Hertz, equivalent to a period of approximately 54 minutes [I]. The structure of Earth's deep 
interior cannot be studied directly. But geologists use seismic waves to determine the depths of layers 
of molten and semi-molten material within Earth. 


P Wave 


Figure 18.18. P-Wave 



The primary wave, or P-wave, is a 
body wave that can propagate 
through the Earth’s core. This 
wave can also travel through 
water. The P-wave is also a 
sound wave. It thus has 
longitudinal motion. Note that 
the P-wave is the fastest of the 
four waveforms. 
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S Wave 



The secondary wave, or S-wave, is a shear wave. It is a type of body wave. The S-wave produces an 
amplitude disturbance that is at right angles to the direction of propagation. Note that water cannot 
withstand a shear force. S-waves thus do not propagate in water. 

Love Wave 



Figure 18.20. Love Wave 

Love waves are shearing horizontal waves. The motion of a Love wave is similar to the motion of a 
secondary wave except that Love wave only travel along the surface of the Earth. Love waves do not 
propagate in water. 
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Rayleigh Wave 



Figure 18.21. Rayleigh Wave, Retrograde 

Rayleigh waves produce retrograde elliptical motion. The ground motion is thus both horizontal and 
vertical. The motion of Rayleigh waves is similar to the motion of ocean waves in 
Figure 18.22 except that ocean waves are prograde. Rayleigh waves resulting from airborne acoustical 
sources may either be prograde or retrograde per Reference [3 I], In some cases, the motion may begin 
as prograde and then switch to retrograde. Airborne acoustical sources include above ground 
explosions and rocket liftoff events. 


Wave Propagation 

-► 



Figure 18.22. Ocean Surface Wave Particle Motion, Prograde 

The Love and Rayleigh waves are both surface waves. These are the two seismic waveforms which can 
cause the most damage to building, bridges and other structures. 

As an aside, seismic and volcanic activity at the ocean floor generates a water-borne longitudinal wave 
called a T-wave, or tertiary wave. These waves propagate in the ocean’s SOFAR channel, which is 
centered on the depth where the cumulative effect of temperature and water pressure combine to 
create the region of minimum sound speed in the water column. SOFAR is short for “Sound Fixing and 
Ranging channel.” These T-waves may be converted to ground-borne P or S-waves when they reach the 
shore. Acoustic waves travel at 1500 m/s in the ocean whereas seismic P and S-waves travel at 
velocities from 2000 to 7000 m/s in the crust. 
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18.6.2 Seismic Response Spectrum Method 

Professors Theodore von Karman and Maurice Biot were very active in the early 1930s in the 
theoretical dynamics aspects of what would later become known as the response spectrum method in 
earthquake engineering. Biot proposed that rather than being concerned with the shape of the input 
time history, engineers should instead use a method describing the response of systems to those shock 
pulses. The emphasis should instead be on the effect, as represented by the response of a series of 
single-degree-of-freedom oscillators, similar to that previously shown for the case of a half-sine input in 
Figure 18.14. 

Practical use of the response spectrum method had to wait until the 1970s due to the intricacy of the 
response calculation for complex, oscillating pulses which required digital computers. Time was also 
needed to establish and publicize databases of strong motion acceleration time histories. 

The response spectrum method was adopted for pyrotechnic shock in the aerospace industry and 
renamed as shock response spectrum. 

18.6.3 El Centro Earthquake 



Figure 18.23. El Centro, Imperial Valley Earthquake Damage 


Nine people were killed by the May 1940 Imperial Valley earthquake. At Imperial, 80 percent of the 
buildings were damaged to some degree. In the business district of Brawley, all structures were 
damaged, and about 50 percent had to be condemned. The shock caused 40 miles of surface faulting on 
the Imperial Fault, part of the San Andreas system in southern California. Total damage has been 
estimated at about $6 million. The magnitude was 7.1. The was the first major earthquake for which 
strong motion acceleration data was obtained that could be used for engineering purposes. 
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North-South Axis 
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Figure 18.24. El Centro Earthquake, Triaxial Time History 
The highest excitation is in the North-South axis, parallel to the ground. 
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Figure 18.25. El Centro Earthquake North-South SRS, Three Damping Cases 

The acceleration levels reached 1.5 G for the 1 % damping curve. Recall that large civil engineering 
structures can have nonlinear damping. The damping values tend to increase as the base excitation 
levels increase as shown for the Transamerica Title Building in Section 10.7. 
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Pseudo Velocity Shock Response Spectrum Q=10 



10' 1 10 ° 10 1 


Natural Frequency (Hz) 

Figure 18.26. El Centro Earthquake Tripartite SRS 

Seismic SRS curves are often plotted in tripartite format which displays relative displacement, pseudo 
velocity and acceleration responses all on the same graph. The pseudo velocity PV is calculated from 
the relative displacement Z as 

PV = <D n Z (18.10) 

Stress can be calculated from pseudo velocity using the methods in Section 19. The acceleration curve 
might be the most important design metric for equipment mounted inside a building. The relative 
displacement might be the most important concern for analyzing the foundational strength. The curves 
also show design tradeoffs. Lowering the building’s natural frequency below 1 Hz reduces the 
acceleration response but increases the relative displacement. Note that Q=10 is the same as 5% 
damping. 
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18.6.4 Golden Gate Bridge 



Figure 18.27. Golden Gate Suspension Bridge, San Francisco, California 


In addition to traffic loading, the Golden Gate Bridge must withstand the following environments: 

1. Earthquakes, primarily originating on the San Andreas and Hayward faults 

2. Winds of up to 70 miles per hour 

3. Strong ocean currents 

The Golden Gate Bridge has performed well in all earthquakes to date, including the 1989 Loma Prieta 
Earthquake. Several phases of seismic retrofitting have been performed since the initial construction. 
The bridge’s fundamental mode is a transverse mode with a natural frequency of 0.055 Hz, with a period 
of 18.2 seconds 

Note that California Department of Transportation (CALTRANS) standards require bridges to 
withstand an equivalent static earthquake force (EQ) of 2.0 G. This level is plausibly derived as a 
conservative envelope of the El Centro SRS curves in Figure 18.25. 
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18.6.5 Vandenberg AFB 



Figure 18.28. Rocket Launch from Vandenberg AFB, California 

The Vandenberg launch site is near the San Andreas fault system. The vehicle is mounted on the pad as 
a tall cantilever beam. The vehicle must be analyzed to verify that it can withstand a major seismic 
event. The vehicle may be mounted to the pad for only two weeks prior to launch. The odds of an 
earthquake occurring to that time window are miniscule. But the launch vehicle and payload together 
may cost well over $1 billion. The risk thus necessitates the analysis. Areas of concern are the loads 
imparted at the launch vehicle’s joints and to the payload. 
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Figure 18.29. NASA SRS Curves for Launch Vehicles at Vandenberg 

SRS curves are given for three damping cases. The curves are taken from Reference [29]. The vehicle 
would typically be analyzed as a multi-degree-of-freedom system via a finite element model. Each SRS 
curve could be applied to the model using a modal combination method. An alternative is to synthesize 
a time history to satisfy a selected SRS curve. The time history could then be applied to the model via a 
modal transient analysis. 


18.6.6 Seismic Testing 



Figure 18.30. Electrical Power Generator Seismic Shock Test 
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The diesel generator is mounted onto a platform at the top of a shaker table which is located below the 
ground floor. This could be an emergency power generator for a hospital in an active seismic zone. A 
video clip of the test is available on YouTube at: https://youtu.be/5uSql7kSYE4 

18.7 Pyrotechnic Shock 


18.7.1 Flight Events 





Figure 18.31. Stage Separation Ground Test, Linear Shaped Charge 

The plasma jet cuts the metal inducing severe mechanical shock energy, but the smoke and fire would 
not occur in the near-vacuum of space. 

Launch vehicle avionics components must be designed and tested to withstand pyrotechnic shock from 
various stage, fairing and payload separation events that are initiated by explosive devices. Solid rocket 
motor ignition is another source of pyrotechnic shock. The source shock energy can reach thousands 
of acceleration Gs at frequencies up to and beyond 100 kHz. The corresponding velocities can reach a 
few hundred in/sec, well above the severity thresholds in the empirical rules-of-thumb in Section 19.4. 

Empirical source shock levels for a variety of devices are given in References [29], [32], [33]. These 
levels are intended only as preliminary estimates for new launch vehicle designs. The estimates should 
be replaced by ground test data once the launch vehicle hardware is built and tested. 
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18.7.2 Pyrotechnic Device Photo Gallery 



Figure 18.32. Metal Clad Linear Shaped Charge 

The chevron focuses a pyrotechnic plasma jet at the launch vehicle’s separation plane material. Severe 
shock levels are generated as a byproduct. 




Figure 18.33. Frangible J oint 

A frangible joint may be used for stage or fairing separation. The key components of a frangible joint: 

• Mild Detonating Fuse (MDF) 

• Explosive confinement tube 

• Separable structural element 

• Initiation manifolds 

• Attachment hardware 

The hot gas pressure from the MDF detonation cause the internal tube to expand and fracture the joint. 
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Figure 18.34. Frangible Nuts, Hold Down Posts 

The purpose of the nuts was to hold the SRBs in place against wind and ground-borne excitation. The 
nuts were separated just before liftoff. 
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Figure 18.35. Clamp Band with Pyrotechnic Bolt-Cutters 
(image courtesy of European space agency) 

Clamp bands are often used for payload separation from launch vehicle adapters. They are also 
commonly used for stage separation in suborbital launch vehicles, similar to the one in Figure 4.1. A 
pyrotechnic bolt-cutter uses an explosive charge to drive a chisel blade to cut the band segments’ 
connecting bolt. The cutters produce some shock energy, but much of the shock is due to the sudden 
release of strain energy in the preloaded clamp band. This action can excite a ring mode in the radial 
axis. Recall Section 7.3. 

The total clamp band release shock is significantly less than linear shaped charge and frangible joint. 
Note that an analysis must be performed to verify that no gapping will occur in the between the band 
and the joint as the vehicle undergoes bending mode vibration during powered flight. 
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18.8 Pyrotechnic Shock Data 

18.8.1 Initial SRS Slopes 



Figure 18.36. Expected Pyrotechnic SRS Initial Slope Limits 

Near-field pyrotechnic shock can be difficult to measure accurately. The accelerometer data may have a 
baseline shift or spurious low frequency transient. This error could be a result of the accelerometer’s 
own natural frequency being excited or to some other problem. Aerospace pyrotechnic shock SRS 
specifications usually begin at 100 Hz due to the difficult in accurately measuring low-amplitude, low- 
frequency shock, while simultaneously measuring high-amplitude, high-frequency shock. 

There are several methods for checking whether the data is acceptable. One is to the check the initial 
slopes of both the positive and negative SRS curves. Each should have an overall trend of 6 to 12 
dB/octave. The actual slopes may have local peaks and dips due to low frequency resonances. A 12 
dB/octave slope represents constant displacement and zero net velocity change. A 6 dB/octave slope 
indicates constant velocity. Recall the slope formulas in Section 16.1.2. 

A second method for checking data accuracy is to verify that the positive and negative SRS curves are 
within about 3 dB of each other across the entire natural frequency domain. A third method is to 
integrate the acceleration time history to velocity. The velocity time history should oscillate in a stable 
manner about the zero baselines. 

These verification goals are challenging to meet with near-field shock measurements of high-energy 
source shock, such as that from linear shaped charge and frangible joints. In practice, some high-pass 
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filtering or spurious trend removal may be necessary. There is no one right way to perform this “data 
surgery.” It is a matter of engineering judgment. 

18.8.2 Re-entry Vehicle Separation, Flight Data 

RV Separation, Flight Accelerometer Data 



Figure 18.37. Re-Entry Vehicle, Separation Shock, Near-Field Measurement, Time History 
The source device was linear shaped charge. 
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Figure 18.38. Re-Entry Vehicle, Separation Shock, Near-Field Measurement, SRS 
The SRS reached 20,040 G at 2400 Hz, which is an extreme level per the rule-of-thumb Section 19.4. 

18.9 Water Impact Shock 



Figure 18.39. Solid Rocket Booster Water Impact Shock 
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Each of the two Space Shuttle Boosters was recovered and refurbished after every flight. Each booster 
contained sensitive avionics components which underwent shock testing for the water impact event. 



Figure 18.40. Solid Rocket Booster Water Impact Shock, Time History 

The data is from the STS-6 Mission. The accelerometer was mounted at the forward end of the booster 
adjacent to a large IEA avionics box. This was the worst-case shock event for this component. 
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Pseudo Velocity Shock Response Spectrum Q=10 



Natural Frequency (Hz) 

Figure 18.41. Solid Rocket Booster Water Impact Shock, Tripartite SRS 


The maximum acceleration response is 257 G at 85 Hz. The maximum pseudo velocity response is 201 
in/sec at 76 Hz, which is severe per the rule-of-thumb Section 19.4. 

18.10 Shock Response Spectrum Synthesis 

18.10.1 Synthesis Objectives 

Consider an SRS specification for a component, subsystem or a large structure. The article should be 
tested, if possible, to verify that it can withstand the shock environment. Certain shock tests can be 
performed on a shaker table, like the generator test in Figure 18.30. This requires synthesizing an 
acceleration time history to satisfy the SRS. The net velocity and net displacement must each be zero 
for this test. These requirements can be met with a certain type of wavelet series, where the individual 
wavelets may be nonorthogonal. The resulting wavelet time history should meet the SRS within 
reasonable tolerance bands, but it may not “resemble” the expected time history which is a limitation of 
this method. Another requirement is that the peak acceleration, velocity and displacement values must 
be within the shaker table’s capabilities. 
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An innovative method for meeting the SRS with a wavelet series that resembles one or more measured 
shock time histories is given in Reference [34]. 

A synthesized time history can also be used for modal transient analysis. This could be done for articles 
which are too large or heavy for shaker tables. This analysis can also be done on small components 
prior to shock testing to determine whether they will pass the test. Or the analysis could be done in 
support of isolator mounting design. There is again a need for the synthesized time history to have net 
velocity and net displacement which are each zero, to maintain numerical stability in stress calculations 
from relative displacement values. 

Damped-sines can be used for modal transient analysis where the goal is to meet the SRS with a time 
history that plausibly resembles the expected field shock event. But damped-sines do not meet the 
desired zero net velocity and displacement goals. The workaround is to first synthesize a damped-sine 
series to meet the SRS and then reconstruct it via a wavelet series, in Rube Goldberg fashion. 

Wavelet and damped-sine synthesis are shown in the following examples. 

18.10.2 Wavelet Synthesis 

18.10.2.1 Wavelet Equation 

A wavelet is a sine function modulated by another sine function. The equation for an individual wavelet 
is 


W m (t) = 


A m sin 
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(18.11) 


where 


W m (t) = acceleration of wavelet mat timet 

A m = wavelet acceleration amplitude 

f m = wavelet frequency 

N m = number of half-sines 

t dm = wavelet time delay 

Note that N m must be an odd integer greater than or equal to 3. This is required so that the net 
velocity and net displacement will each be zero. 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


225 












enDAQ 


Chapter: Shock Testing & Analysis 


The total acceleration at time t for a set of n wavelets is 



Zw,( t ) 

m=l 


( 18 . 12 ) 


WAVELET 1 FREQ = 74.6 Hz 
NUMBER OF HALF-SINES = 9 DELAY = 0.012 SEC 



Figure 18.42. Sample Wavelet 

A sample, individual wavelet is shown in the above figure. This wavelet was a component of a previous 
analysis for an aerospace project. 

18.10.2.2 Wavelet Synthesis Example 

Consider the specification: MIL-STD-810E, Method 516.4, Crash Hazard for Ground Equipment. 

Table 18.3. SRSQ=10, Crash Hazard Specification 


Natural 

Frequency (Hz) 

Peak Accel 

(G) 

10 

9.4 

80 

75 

2000 

75 
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Synthesize a series of wavelets as a base input time history for a shaker shock test to meet the Crash 
Hazard SRS. The goals are: 

1. Satisfy the SRS specification 

2. Minimize the displacement, velocity and acceleration of the base input 
The synthesis steps are shown in the following table. 


Table 18.4. Wavelet Synthesis Steps 


Step 

Description 

i 

Generate a random amplitude, delay, and half-sine number for each wavelet. 

Constrain the half-sine number to be odd. These parameters form a wavelet table. 

2 

Synthesize an acceleration time history from the wavelet table. 

3 

Calculate the shock response spectrum of the synthesis. 

4 

Compare the shock response spectrum of the synthesis to the specification. Form a 
scale factor for each frequency. 

5 

Scale the wavelet amplitudes. 

6 

Generate a revised acceleration time history. 

7 

Repeat steps 3 through 6 until the SRS error is minimized or an iteration limit is 
reached. 

8 

Calculate the final shock response spectrum error. 

Also calculate the peak acceleration values. 

Integrate the signal to obtain velocity, and then again to obtain displacement. 

Calculate the peak velocity and displacement values. 

9 

Repeat steps 1 through 8 many times. 

10 

Choose the waveform which gives the lowest combination of SRS error, 
acceleration, velocity and displacement. 


The resulting time history and SRS are shown in Figure 18.43 and Figure 18.44, respectively. 
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Acceleration 



Velocity 




Time (sec) 


Figure 18.43. Crash Flazard Time History Synthesis 

The acceleration time history has a reverse sine sweep character. It is an efficient and optimized 
waveform for a shaker shock test, and it satisfies the SRS as shown in the next figure. A drawback is 
that it does not resemble an actual crash shock time history. 
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Figure 18.44. Crash Hazard SRS 

The positive and negative curves are from the synthesized waveform. The tolerance bands are set at +3 dB. 


18.10.3 Damped-Sine Synthesis 

18.10.3.1 Damped-Sine Equation 

The equation for an individual damped sinusoid is 


w n (t) = 


0, for t < t 


dn 


A n expr-^ n 0 n (t-t dn )lsinrro n (t-t dn )l, for t >t dn 


where 


W n (t) = acceleration of damped sinusoid nat timet 

A n = acceleration amplitude 

= angular frequency 
= damping value 
tdn = time delay 


(18.13) 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


229 



















enDAQ 


Chapter: Shock Testing & Analysis 


The total acceleration at time t for a set of damped sinusoids is 

N 

x(t) = Z W n(t) (18.14) 

n=l 


18.10.3.2 Damped-Sine Example 

Consider the following specification which could represent a stage separation shock level as some 
location in a launch vehicle. A modal transient finite element analysis is to be performed on a 
component to verify that the component will pass its eventual shock test. The immediate task is to 
synthesize a time history to satisfy the SRS. The time history should “resemble” an actual pyrotechnic 
shock pulse. 

Note that pyrotechnic SRS specifications typically begin at 100 Hz. The author’s rule-of-thumb is to 
extrapolate the specification down to 10 Hz in case there are any component modes between 10 and 
100 Hz. This guideline is also to approximate the actual shock event which should have an initial ramp 
somewhere between 6 and 12 dB/octave. 

Table 18.5. SRSQ=10, Stage Separation Shock 


Natural 

Frequency (Hz) 

Peak Accel 

(G) 

10 

10 

2000 

2000 

10,000 

2000 


The specification has an initial slope of 6 dB/octave. The synthesis steps are shown in the following table. 
Table 18.6. Damped-Sine Synthesis Steps 


Step 

Description 

i 

Generate random values for the following for each damped sinusoid: amplitude, 
damping ratio and delay. The natural frequencies are taken in one-twelfth octave 
steps. 

2 

Synthesize an acceleration time history from the randomly generated parameters. 

3 

Calculate the shock response spectrum of the synthesis. 

4 

Compare the shock response spectrum of the synthesis to the specification. Form a 
scale factor for each frequency. 

5 

Scale the amplitudes of the damped sine components. 

6 

Generate a revised acceleration time history. 
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7 

Repeat steps 3 through 6 as the inner loop until the SRS error diverges. 

8 

Repeat steps 1 through 7 as the outer loop until an iteration limit is reached. 

9 

Choose the waveform which meets the specified SRS with the least error. 

10 

Perform wavelet reconstruction of the acceleration time history so that velocity and 
displacement will each have net values of zero. 


The resulting time history and SRS are shown in Figure 18.45 and Figure 18.46, respectively. 





Figure 18.45. Stage Separation Time History Synthesis 
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The acceleration time history somewhat resembles a mid or far-field pyrotechnic shock. The velocity 
and displacement time histories each have a stable oscillation about their respective baselines. 



Figure 18.46. Stage Separation SRS 
The positive and negative curves are from the damped-sine synthesis. 


18.10.4 Modal Transient Finite Element Analysis for Uniform Base Excitation 

Consider a rectangular plate mounted to a base at each of its four corners. The plate is to be subject to 
uniform seismic excitation. There are two methods to apply the base excitation in a finite element 
analysis, as shown in Figure 18.47 and Figure 18.48. 



u 1 = t^=u 3 =u 4 


Figure 18.47. Direct Enforce Acceleration Method 
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The direct enforcement method is computationally intensive, requiring matrix transformations and a 
matrix inversion as shown in Reference [35]. 



Figure 18.48. Seismic Mass Method 

The seismic mass is chosen to be several orders of magnitude higher in mass than the plate. An 
equivalent force is calculated and apply to the seismic mass to excite the desired acceleration at each of 
the plate’s corners. This method adds a degree-of-freedom to the plate system resulting in a rigid-body 
mode at zero frequency. But the remaining natural frequencies and mode shapes should be the same as 
if the plate were mounted normally to its joining structure. The author’s experience is that the seismic 
mass method is faster and more accurate and reliable than the direct enforced method. This method 
was introduced in Section 12.3.3. 

18.10.5 Shock Fields 

18.10.5.1 Near-Field 

The near-field environment is dominated by direct stress wave propagation from the source. Peak 
accelerations in excess of 5000 G occur in the time domain with a frequency content extending beyond 
100 kHz. The near-field usually includes structural locations within approximately 15 cm of the source 
for severe devices such as linear shaped charge and frangible joint. No shock-sensitive hardware should 
be mounted where it would be exposed to a near-field environment. 

18.10.5.2 Mid-Field 

The mid-field environment is characterized by a combination of wave propagation and structural 
resonances. Peak accelerations may occur between 1000 and 5000 G, with substantial spectral content 
above 10 kHz. The mid-field for intense sources usually includes structural locations between 
approximately 15 and 60 cm of the source, unless there are intervening structural discontinuities. 
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18.10.5.3 Far-Field 

The far-field environment is dominated by structural resonances. The peak accelerations tend to fall 
below 2000 G, and most of the spectral content below 10 kHz. The far-field distances occur outside the 
mid-field. The typical far-field SRS has a knee frequency corresponding to the dominant frequency 
response. The knee frequency is the frequency at which the initial ramp reaches the plateau in the log- 
log SRS plot. 


18.10.5.4 Joint & Distance Attenuation 

The source shock energy is attenuated by intervening material and joints as it propagates from the near¬ 
field to the far-field. Empirical distance and joint attenuation factors for the SRS reduction are given in 
References [29], [32], [33]. A typical attenuation curve from [32] is given in Figure 18.49, which 
assumes an input source shock SRS consisting of a ramp and plateau in log-log format. Such curves 
should be used with caution given that the attenuation is highly dependent on damping and structural 
details. The curves can be used for preliminary estimates for new launch vehicle designs, but ground 
separation tests are still needed for the actual launch vehicle hardware. These tests are needed to 
measure the source shock as well as the levels at key component mounting locations. 



Figure 18.49. Shock Response Spectrum Versus Distance from Pyrotechnic Shock Source 
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18.10.6 Shock Mitigation 



Figure 18.50. Sensor Electronics, Wire Rope Isolators 
(image courtesy of NASA/J PL) 

The source shock energy is attenuated as it propagates to avionics mounting locations through the 
launch vehicle’s material and joints. The input shock to a component can be mitigated by mounting the 
component as far away from the source device as possible. Another effective attenuation technique is 
to mount the component via elastomeric bushings or wire rope isolators. The NASA Mars Science 
Laboratory Sensor Support Electronics unit is mounted on vibration isolators in Figure 18.50. 



isolator Bustling 


Wooden Avionics Shell 


Figure 18.51. SCUD-B Avionics Component Isolation 
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The bushings are made from some type of rubber or elastomeric compound. The bushings provide 
damping, but their main benefit is to lower the natural frequency of the system. The isolators thus 
attenuate the shock and vibration energy which flows from the instrument shelf into the avionics 
component. 

18.10.7 Pyrotechnic Shock Testing Methods 



Figure 18.52. Near-Field Shock Simulation using a Plate 

The test component is mounted on other side of plate. The source device is a textile explosive cord 
with a core load of 50 gr/ft (PETN explosive). Up to 50 ft of Detonating Cord has been used for some 
high G tests. The maximum frequency of shock energy is unknown so analog anti-aliasing filters are 
needed for the accelerometer measurements per the guidelines in Section 14. 

Note that a component may be mounted in the mid-field shock region of a launch vehicle. The SRS test 
level derivation for this zone may include a significant dB factor for uncertainty or as a qualification 
margin. This conservatism may require a near-field-type shock test for a component that is actually 
located in a mid-field zone. This situation can also occur for components mounted in far-fields. 
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Figure 18.53. NASA J PL Tunable Shock Beam 

The NASA/JPL Environmental Test Laboratory developed and built a tunable beam shock test bench 
based on a design from Sandia National Laboratory many years ago. The excitation is provided by a 
projectile driven by gas pressure. The beam is used to achieve SRS specifications, typically consisting of 
a ramp and a plateau in log-log format. The intersection between these two lines is referred to as the 
“knee frequency.” The beam span can be varied to meet a given knee frequency. The high frequency 
shock response is controlled by damping material. 
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18.10.8 Shock Failure Modes 



Figure 18.54. Sensitive Electronic Parts 

Pyrotechnic shock can cause crystal oscillators to shatter. Large components such as DC-DC 
converters can detached from circuit boards. In addition, mechanical relays can experience chatter or 
transfer. 



Figure 18.55. Shock Test Case History, Adhesive & Solder J oint Failure 

The image shows adhesive failure and rupture of solder joint after a stringent shock test. A large 
deflection of the PCB resulting from an insufficient support/reinforcement of the PCB combined with 
high shock loads can lead to these failures. Staking is needed for parts weighing more than 5 grams. 
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Figure 18.56. Shock Test Case History, Lead Wire Failure 

The image shows a sheared lead between solder joint and winding of coil. The lacing cord was 
insufficient by itself. The lacing should be augmented by staking with an adhesive. 
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19 Pseudo Velocity & the Stress-Velocity Relationship 

19.1 Pseudo Velocity 

The Shock Response Spectrum (SRS) models the peak response of a single-degree-of-freedom (SDOF) 
system to a base acceleration, where the system's natural frequency is an independent variable. The 
absolute acceleration and the relative displacement of the SDOF system can be readily calculated. The 
relative velocity can also be calculated, but the common practice is to approximate the relative velocity 
in terms of the "pseudo velocity." 

The preferred method for calculating the peak pseudo velocity is to multiply the peak relative 
displacement by the angular natural frequency as shown in equation (18.10). An alternate method is to 
divide the peak acceleration by this frequency. 

19.2 Stress-Velocity Relationship 

Howard A. Gaberson (1931-2013) championed the idea that 
dynamic stress correlated better with pseudo velocity than with 
either acceleration or relative displacement. Gaberson was a 
shock and vibration specialist who was with the U.S. Navy Civil 
Engineering Laboratory and later the Facilities Engineering Service 
Center from 1968 to 2000, mostly conducting dynamics research. 
He specialized in shock and vibration signal analysis and published 
more than 100 papers and articles. Gaberson pointed out that the 
kinetic energy in a mechanical system is proportional to velocity 
squared. The pseudo velocity is thus a measure of the stored 
peak energy in the system at a particular frequency and has a 
direct relationship to the survival or failure of the system. The 
following equations are taken from Reference [36]: 

Consider the stress-velocity relationship for a traveling wave in an infinite rod. 

Compression zone —^ Rarefaction zone —^ 



Direction of travel 


Figure 19.1. Propagating Wave in an Infinite Rod 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


240 










enDAQ 


Chapter: Pseudo Velocity & the Stress-Velocity Relationship 


The stress a is proportional to the velocity as 

a(x,t) = —pcv(x,t) ( 19 . 1 ) 

The velocity depends on natural frequency, but the stress-velocity relationship does not. The lack of a 
frequency term in equation (19.1) was another reason that Gaberson promoted the stress-velocity 
relationship. 



Figure 19.2. Propagating Wave in a Finite Rod 

The principle can be extended to a finite rod with common boundary conditions where the peak modal 

stress K] max is 

[ a n ] max — P C [^ n ] max (19.2) 

The index n represents the mode number. The peak stress is mainly due to the fundamental mode in 
most cases. But higher modes can contribute significant stress in certain cases, such as pyrotechnic 
shock response. 

19.3 Beam Stress 

Consider the bending vibration of a simply-supported beam. The modal stress due to bending is 
proportional to the modal velocity as follows. 

[ a n] m ax = c V EA P/ r [ V n] max (19.3) 

The distance from the beam’s neutral axis to the outer fiber is c . 

Note this equation applies to other boundary condition cases per Reference [37]. 

Equation (19.3) can be simplified as follows: 

[ a n ] max — ^P c [ Yn ] max (19.4) 
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Values for the k constant for typical cross-sections are: 


Cross-section 

k 

Solid Circular 

2 

Rectangular 



The method can also be used for complex structures, but with some uncertainty. The stress-velocity 
formula is 

[ a n]max — Kpc [V n ] max ( 19 . 5 ) 

The constant of proportionality K depends on geometry, etc. Bateman [38] claimed the following range 
for complex equipment 4 < K < 8. Gaberson [39] gave a larger range with 1 < K < 10. 

Equation (19.5) could thus be used to “bound the problem” by assuming a high and low K value. 

Note that the stress-velocity equations do not require that peak stress and peak velocity occur at the 
same location. Any stress concentration factor must be applied separately. The equations are useful 
as quick estimates prior to any advanced modeling such as finite element analysis. The equations can 
also be used to “sanity check” any finite element results. 

19.4 Shock Velocity Threshold Rule-of-Thumb 

An empirical rule-of-thumb in MIL-STD-810E [26] and SMC-TR-06-11 [40] states that a shock response 
spectrum is considered severe only if one of its components exceeds the level 

Threshold =[0.8 (G/ Hz) * Natural Frequency (Hz) ] (19.6) 

For example, the severity threshold at 100 Hz would be 80 G. 

This rule is effectively a velocity criterion. MIL-STD-810E states that it is based on unpublished 
observations that military-quality equipment does not tend to exhibit shock failures below a shock 
response spectrum velocity of 100 inches/sec (254 cm/sec). 

Equation (19.6) actually corresponds to 50 inches/sec. It thus has a built-in 6 dB margin of conservatism. 

The 100 inches/sec threshold is defined in part by the observation that the severe velocities which cause 
yield point stresses in mild steel beams turn out to be about 130 inches/sec. 
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19.5 Velocity Threshold from Yield Stress 

The velocity threshold V th can be calculated from the yields stress as follow for a beam with rectangular 
cross section. 


V th = 


^yield 

Kpc 


The parameters for mild steel are given in the following table. 

Table 19.1. Mild Steel Properties 


^ yield 

33 ksi 

p 

0.00075 lbf sec^2/ in^4 

C 

199,000 inches/sec 


(19.7) 


V th = 


33,000 lbf / in 2 


y/3 (0.00075 lbf sec 2 / in 4 1(199,000 in/sec) 


= 128 in/sec 


(19.8) 


Gaberson [39] gave the limits in the following table. 


Table 19.2. Severe Velocities, Fundamental Limits to Modal Velocities in Structures 


Material 

E (psi) 

^yield 

(psi) 

P 

(lbm/in 3 ) 

Rod 

Vmax 

(ips) 

Beam 

Vmax 

(ips) 

Plate 

Vmax 

(ips) 

Douglas Fir 

1.92e+06 

6450 

0.021 

633 

366 

316 

Aluminum 

6061-T6 

10.0e+06 

35,000 

0.098 

695 

402 

347 

Magnesium 

AZ80A-T5 

6.5e+06 

38,000 

0.065 

1015 

586 

507 

Structural Steel, 
High Strength 

29e+06 

33,000 

100,000 

0.283 

226 

685 

130 

394 

113 

342 
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19.6 Morse Chart 



Natural Frequency (Hz) 

Figure 19.3. Morse Chart, Shock Response Spectra, Damage Thresholds 

The Morse Chart in the above figure was derived for aerospace components subjected to shock 
environments [41]. The curves are defined by the following formulas. 


Table 19.3. Morse Chart Severity Thresholds 


Threshold 

Formula 

300 ips 

[ 4.8 (G/Hz) * Natural Frequency (Hz) ] 

100 ips 

[ 1.6 (G/Hz) * Natural Frequency (Hz) ] 

50 ips 

[ 0.8 (G/Hz) * Natural Frequency (Hz) ] 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


244 













enDAQ 


Chapter: Digital Filtering 


20 Digital Filtering 

20.1 Filtering Introduction 

Filtering is a tool for resolving signals. It can be performed on either analog or digital signals. Analog 
anti-aliasing filtering was covered in Section 14.4. The present section is limited to digital filtering with 
an emphasis on the Butterworth filter, which is implemented in the time domain as a digital recursive 
filtering relationship. 

20.2 High-pass & Low-pass Filters 

A high-pass filter is a filter which allows the high-frequency energy to pass through. It is thus used to 
remove low-frequency energy from a signal. 

A low-pass filter is a filter which allows the low-frequency energy to pass through. It is thus used to 
remove high-frequency energy from a signal. 

A band-pass filter may be constructed by using a high-pass filter and low-pass filter in series 

20.3 Butterworth Filter Characteristics 


BUTTERWORTH LOWPASS FILTER L=ORDER Q, = 1 rad/sec 



Figure 20.1. Butterworth Filter, Three Orders, Transfer Function Magnitudes 

A Butterworth filter is one of several common infinite impulse response (MR) filters. Other filters in this 
group include Bessel and Chebyshev filters. These filters are classified as feedback filters. The 
Butterworth filter can be used either for high-pass, low-pass, or band-pass filtering. It is characterized 
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by its cut-off frequency. The cut-off frequency is the frequency at which the corresponding transfer 
function magnitude is -3 dB, equivalent to 0.707. The transfer function curves in Figure 20.1 each pass 
through the same cut-off frequency point. 

A Butterworth filter is also characterized by its order, with the sixth-order preferred in this document. 
A property of Butterworth filters is that the transfer magnitude is -3 dB at the cut-off frequency 
regardless of the order. Other filter types, such as Bessel, do not share this characteristic. 

Consider a lowpass, sixth-order Butterworth filter with a cut-off frequency of 100 Hz. The 
corresponding transfer function magnitude is given in Figure 20.2. 




Figure 20.2. Butterworth Low-Pass Filter, 100 Hz, Transfer Function 

Note that the curve in the above figure has a gradual roll-off beginning at about 70 Hz. Ideally, the 
transfer function would have a rectangular shape, with a corner at (100 Hz, 1.0). This ideal is never 
realized in practice due to stability concerns in the time domain. Thus, a compromise is usually required 
to select the cut-off frequency. The transfer function could also be represented in terms of a complex 
function, with real and imaginary components. A transfer function magnitude plot for a sixth-order 
Butterworth filter with a cut-off frequency of 100 Hz as shown in Figure 20.3. 
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Figure 20.3. Butterworth High-Pass Filter, 100 Hz, Transfer Function 
The curves in the previous figure suggests that filtering could be achieved as follows: 

1. Take the Fourier transform of the input time history 

2. Multiply the Fourier transform by the filter transfer function, in complex form 

3. Take the inverse Fourier transform of the product 

The above frequency domain method is valid. Nevertheless, the filtering algorithm is usually 
implemented in the time domain for computational efficiency, maintaining stability, avoid leakage error, 
etc. 
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Input X 


Time domain 
equivalent of 

H(co) 


Response y 


Figure 20.4. Filter Block Diagram 


The filter transfer function can be represented by H(eo) as shown in Figure 20.4. 

The filtering transfer function is implemented as a digital recursive filtering relationship in the time 
domain. The response y k at time index k is 


Yk = 


Z b n x k-n 

n=0 


L 

Z a nYk-n 


n=0 


( 20 . 1 ) 


The filter order is L. The filter coefficients a n and b n can be derived using the method is Reference 
[42]. The equation is recursive because the output at any time depends on the output at previous times. 

20.4 Phase Correction 



Figure 20.5. Refiltering for Phase Correction 


An ideal filter should provide linear phase response. This is particularly desirable if shock response 
spectra calculations are required. Butterworth filters, however, do not have a linear phase response. 
Other MR filters share this problem. Methods are available, however, to correct the phase response. 
One is the refiltering method in Reference [42]. An important note about refiltering is that it reduces 
the transfer function magnitude at the cut-off frequency to -6 dB, as shown in example in Figure 20.6. 
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180 


Butterworth Lowpass Filter 6th order Refiltering fc=100 Hz 




Figure 20.6. Butterworth Filter 6 th Order Low-Pass Filter with Refiltering, Transfer Function 

20.5 PSD Calculation via Digital Filter 

Recall the synthesized time history from Figure 16.10 which satisfied the Navmat P-9492 PSD 
specification. Calculate the power spectral density by applying band-pass filtering over successive bands. 
This method is useful as a means of better understanding the PSD function. The filtering is performed 
using a Butterworth sixth-order filter without refiltering for phase correction. The bands are full 
octave. The results are shown in the following table. 
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Table 20.1. Band-pass Filter Results, Octave Bands 


fl (Hz) 

fc (Hz) 

fu (Hz) 

Af (Hz) 

GRMS 

GRMS A 2 

GRMS A 2/Hz 

14 

20 

28 

14 

0.36 

0.13 

0.0091 

28 

40 

57 

29 

0.80 

0.64 

0.0226 

57 

80 

113 

56 

1.47 

2.17 

0.0384 

113 

160 

226 

113 

2.16 

4.67 

0.0413 

226 

320 

453 

227 

2.92 

8.54 

0.0377 

453 

640 

905 

452 

3.11 

9.70 

0.0214 

905 

1280 

1810 

905 

3.02 

9.15 

0.0101 

1810 

2560 

3620 

1810 

1.31 

1.71 

0.0009 


The band center frequency is fc. The lower and upper band limits are fl and fc, respectively. The 
bandwidth is Af which is the difference between the band limits. A GRMS value is calculated for each 
band. Then the GRMS value is squared and divided by the bandwidth. The last column is the PSD in 
GRMS A 2/Hz. The time histories for the two rows in the table are shown in Figure 20.7. The 
GRMS A 2/Hz points are plotted along with Navmat PSD in Figure 20.8. Recall that GRMS A 2/Hz is 
abbreviated as G A 2/Hz. The bands are in full octave format. 
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Figure 20.7. Synthesized Time History with Two Filtering Cases, Close-up View 


BP Filtered, 160 Hz Band 


i-1-1-1-1-1-1-1-r 



J_I_I_I_I_I_I_I_L 


Synthesis 



3 5.1 5.2 5.3 5.4 5.5 5.6 5.7 5.8 5.9 6 

BP Filtered, 40 Hz Band 

- 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 




_i_i_i_i_i_i_i_i_i_ 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


251 

























enDAQ 


Chapter: Digital Filtering 
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Frequency (Hz) 

Figure 20.8. PSD Comparison 

The PSD points calculated via band-pass filtering track the specification fairly well. The dropout for the 
last point is not a concern because the bandwidth extended from 1810 to 3620 Hz, but the specification 
stopped at 2000 Hz. 


20.6 Seismic Data Filtering 

The following example shows how filtering can aid in analyzing seismic data from a seismometer. 
Filtering is used to find onset of P-wave in seismic time history from Solomon Island earthquake, 
Magnitude 6.8, October 8, 2004. The measured data is from a homemade seismometer in Mesa, 
Arizona. 
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Figure 20.9. Author's Homemade Lehman Seismometer 

The boom is a horizontal pendulum. It has a period of 14.2 seconds, equivalent to a natural frequency of 
0.071 Hz. A sensor at the free end measures the displacement. The boom length is 64 inch. The total 
frame height is 35 inch. The boom has a knife edge that pivots against a bolt head in the lower cross¬ 
beam of the frame. 

The boom is suspended from the frame by a wire cable. The cable is attached to the top cross-beam of 
the frame. The other end of the cable is attached to the boom, about two-thirds of the distance from 
the pivot to the free end of the boom. The pivot point is offset from the top cable attachment point. 
Thus, the boom oscillates as if it were a "swinging gate." 

The plate supporting the frame has three adjustable mounting feet. The feet can be adjusted to tune the 
pendulum to the desired natural frequency. Furthermore, the wire cable has a turnbuckle which is used 
to adjust height of the free end of the boom. The detached frame in the center of the figure is used for 
assembly and to limit the displacement during tuning. 
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Figure 20.10. Non-contact Relative Displacement Sensor 

The classic sensor for Lehman seismometers is a magnetic coil attached to the boom. The coil moves 
between the poles of a magnet, thus inducing a voltage proportional to the velocity. The signal requires 
tremendous amplification. The author has better success with the inductive displacement sensor. 



Figure 20.11. Ballast Mass Partially Submerged in Oil 
The damping method is shown in this figure. 
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Figure 20.12. Pivot at the End of the Boom 
The chisel blade butted up against a chrome-plated bolt head. 


Seismometer Calibration T=14.2 sec, 9.8% Damping 



-5 - 1 - 1 - 1 - 1 - 
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Figure 20.13. Seismometer Free Vibration Response 

The seismometer was given an initial displacement and then allowed to vibrate freely. The period was 
14.2 seconds, with 9.8% damping. The corresponding frequency is 0.0705 Hz. 
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LEHMAN SEISMOMETER HORIZONTAL RESPONSE TO 
SOLOMON ISLAND EARTHQUAKE UTC 2004/10/08 03:27:53 



TIME (SEC) 

Figure 20.14. Measured Solomon Earthquake Time History, Lehman Seismometer 

The trace in the above figure is the displacement time history of the Solomon Islands earthquake on 
October 8, 2004 as measured by the Lehman seismometer via the inductive sensor. The seismometer 
was located at home in Mesa, Arizona. The data was acquired by a Nicolet Vison system. The Nicolet 
sample rate was set to 50 samples per second with its lowpass filter set to 5 Hz. 

A Krohn-Hite filter, model 3343, was used to high-pass filter the analog displacement signal at 0.03 Hz 
prior to its input to the Nicolet system. It also provided a 20 dB gain. 

The time is referenced to the earthquake occurrence using the USGS data. The plot’s Y-axis is labeled 
as relative displacement because it is the response of the boom relative to the ground. Further 
calculation would be required to estimate the true ground motion. The time history shows that the 
Earth is remarkably reverberant. The oscillations last well over one hour. The phase components are P 
primary wave, S secondary or shear wave, and LQ Love wave. Recall the seismic wave diagrams in 
Section 18.6.1. 

The P-wave is indiscernible against the background microseismic noise. Nevertheless, it can be extracted 
by additional filtering, as shown in the next figure. 
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Figure 20.15. High-pass Filtered Seismic Waveform Showing Onset of P-wave 


The P-wave arrives at in Mesa, Arizona about 760 seconds after it was generated near the Solomon 
Islands. The S-wave arrives around 1400 seconds. There are also some intermediate waveforms which 
begin as P-waves and transform into other types through reflection and refraction. Both the P and S- 
waves are body waves which travel faster than surface waves. The frequency content of body and 
surface waves is shown in the following table, as taken from Reference [43]. 


Table 20.2. Seismic Waveform Periods and Frequencies 


Wave Type 

Period (sec) 

Natural Frequency (Hz) 

Body 

0.01 to 50 

0.02 to 100 

Surface 

10 to 350 

0.003 to 0.1 


High-Pass Filtered Time History 0.2 Hz 



20.7 Integration and Differentiation of Accelerometer Signals 

There are occasional needs to integrate accelerometer time histories to velocity. This may be done to 
evaluate the accuracy of pyrotechnic shock data where the velocity is expected to have a stable 
oscillation about its zero baseline. Accelerometer data may also be double-integrated to displacement, 
but the accuracy depends on the frequency response of the accelerometers. Some accelerometer types 
such as variable capacitance and servo motor designs can accurately measure acceleration down to zero 
frequency. The common piezoelectric accelerometer may have a practical lower limit of a few Hertz, 
however. Some type of high-pass filter is needed in this case. The high-pass filtering may be performed 
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using an analog filter in the accelerometer’s signal conditioner. Or the filtering may be performed 
digitally as a post-processing step. 

Another integration scenario is the case where an acceleration time history is synthesized for use in a 
modal transient analysis. The corresponding velocity and displacements should each be stable, especially 
if dynamic stresses are to be calculated. Consider the simple case of the synthesized white noise 
acceleration time history in Figure 20.16. 

The integration steps are performed using the trapezoidal method, which has sufficient accuracy. The 
resulting velocity displacement each has a drift. These effects arise in part because the integration does 
not account for initial velocity and initial displacement, each of which are either unknown or unspecified. 
The next task is to correct the acceleration time history so that its velocity and displacement will each 
be stable. There is no one right way to do this, but an example is shown in Figure 20.17. The 
correction steps were: 

1. Remove mean from acceleration 

2. High-pass filter acceleration at 6 Hz, Butterworth 6th order with refiltering for phase correction 

3. Apply fade in and out to acceleration using 2% of the total duration at each end 

4. Integrate to velocity 

5. Remove mean from velocity 

6. Apply fade in and out to velocity using 2% of the total duration at each end 

7. Integrate to displacement 

8. Perform first-order trend removal on displacement 

9. Apply fade in and out to displacement using 2% of the total duration at each end 

10. Differentiate displacement to velocity 

I I. Differentiate velocity to acceleration Integrate corrected acceleration to velocity 

12. Integrate velocity to displacement 

13. Verify that the resulting velocity and displacement each have stable oscillations about their 
respective zero baselines and that each begins and ends at zero 
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Figure 20.16. White Noise Acceleration Time History with Corresponding Velocity and Displacement 
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Figure 20.17. Acceleration Time History Correction and Verification 

Nine subplots appear as follows: 

Row 1, Col 1 - processed acceleration prior to integration 

Row 2, Col 1 - integrated velocity 

Row 3, Col 1 - double integrated displacement 

Row 3, Col 2 - double integrated displacement (repeated) 

Row 2, Col 2 - differentiated velocity from double integrated displacement 

Row 1, Col 2 - double differentiated acceleration from double integrated displacement 
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Row 1, Col 3 - double differentiated acceleration (repeated) 

Row 2, Col 3 - integrated velocity from double differentiated acceleration 

Row 3, Col 3 - double integrated displacement from double differentiated acceleration 

The corrected acceleration is the time history shown in Row 1, Col 2 and repeated in Row 1, Col 3. 



Figure 20.18. Fade In and Out Windows 

Note that the fade in and out windows are similar to those shown in the above figure. 

20.8 Band Limited White Noise 


White Noise 



Figure 20.19. White Noise Segment Before and After Low-Pass Filtering 
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White noise is an idealized concept. True white noise would have a flat PSD level from zero to infinity 
frequency units. Digital white noise is inherently frequency limited by its sample rate. Furthermore, 
white noise may be subject to low-pass filtering to smooth it prior to its use as an excitation function, as 
shown in Figure 20.19. The top plot is a white noise sampled at 10,000 samples/sec and has a choppy 
appearance. Aliasing would be suspected if the top plot were measured data. 

The bottom plot is the same time history low-pass filtered at 2000 Hz which smooths the data. Either 
time history should be suitable for either applied force or base excitation of an SDOF system with a 
natural frequency of 200 Hz, for example. Recall that the SDOF system itself behaves a mechanical filter. 
But the low-pass filtered version could have some potential numerical accuracy benefits, although an 
investigation of this hypothesis is left as future work. 


20.9 


Launch Vehicle 


Wind Velocity 



Coupled-Loads Analysis Filtering 



Figure 20.20. Sources of Aero-Mechanical Loads 

Coupled Loads Analysis (CLA) predicts payload & launch vehicle responses due to major dynamic and 
quasi-static loading events. This analysis is performed prior to launch. It can also be performed as post¬ 
flight data reconstruction using flight accelerometer data. Flight accelerometer data is low-pass filtered 
for coupled-loads analyses. The cut-off frequency varies by launch vehicle, payload, key events, etc. The 
primary sources of these low frequency loads are: 

1. Pre-launch events: ground winds, seismic loads 

2. Liftoff: engine/motor thrust build-up, ignition overpressure, pad release 

3. Air-loads: buffet, gust, static-elastic 

4. Liquid engine ignitions and shutdowns 

Reference [44] gives the following frequency guidelines for CLA. The low-frequency dynamic response 
of the launch vehicle/payload system to transient flight events is typically analyzed from 0 Hz to 100 Hz. 
The response frequency domain can be up to 150 Hz for some small launch vehicles. 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


262 






enDAQ 


Chapter: The Great Spacecraft Base Input Vibration Test Debate 


21 The Great Spacecraft Base Input Vibration Test Debate 

21.1 Debate Introduction 



Figure 21.1. Shaker Table Vibration Test 


Engineers prepare the MESSENGER spacecraft for a vibration test on a shaker table at The Johns 
Hopkins University Applied Physics Laboratory, Laurel, Maryland. The purpose of this test is to verify 
that the spacecraft can withstand the launch vehicle environment during liftoff and ascent. 

Flight configured spacecraft are subjected to base input vibration tests for certain programs. The 
spacecraft are often one-of-a-kind, so the vibration test is effectively a proto-qualification test covering 
both design and workmanship verification. The tests may be sinusoidal or random. The sine vibration is 
typically at low frequencies, below 100 Hz. The author’s colleagues are divided on whether these 
spacecraft system-level tests are prudent and effective. 
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21.2 Arguments against Spacecraft Vibration Testing 

The following assertions are made by A.M. Kabe and E. Perl from The Aerospace Corporation in 
informal working groups meetings. 

1. Vibration tables cannot replicate the impedance of the launch vehicle interface, nor the 
interaction that occurs between the launch vehicle and spacecraft when they are a coupled 
system; hence, the modes of vibration will not be the same as in flight. 

2. Only translational motions are applied at the base, one axis at a time, whereas during flight, 
the launch vehicle/spacecraft system will vibrate simultaneously in all six degrees of 
freedom at each mass point and at each interface point between the launch vehicle 
and spacecraft. 

3. The total acceleration load during powered flight also depends on the spacecraft rigid-body 
acceleration which a shaker cannot replicate. 

4. Derivation of a “base input” environment from a few accelerometer locations at the 
launch vehicle/spacecraft interface will generally lead to an over prediction of the motions 
at the interface, since local deformations are mapped on the assumption that the interface 
acts as a rigid plane. 

5. The use of (response)/(base motion) ratios to extract damping, a common practice, is not a 
valid approach for multi-degree-of-freedom systems it fails to account for the mode 
participation factor. 

6. The test article may not include the actual spacecraft launch vehicle adapter or the 
propellant mass in the tanks because of safety and contamination concerns. 

7. The test requirement forces the spacecraft organization to design its system to not only 
survive the launch environment, but also to survive an artificial test that more often than 
not produces overly conservative loads in many parts of the structure while not adequately 
testing others. 

8. The test can pose unnecessary risk of damaging flight hardware late in the program. 


Note that A.M. Kabe advocates acoustic reverberant chamber testing of spacecraft as a workmanship 
screen, as an alternative to base shake testing. He also favors shaker table vibration testing on a 
component or subsystem level. 

21.3 Arguments for Spacecraft Vibration Testing 

The following justification points are made by NASA engineers Daniel Kaufman, Scott Gordon, Steve 
Hendricks and Dennis Kern. 

1. Essentially all current launch vehicle organizations (Delta, Atlas, Taurus, Pegasus, Ariane, 

Hll, Proton, Long March, Falcon, etc.) specify and require or strongly recommend a 
spacecraft sine or random vibration test. 

2. Testing is also required by NASA documents, such as References [45] and [46]. 

3. Insurance companies require vibration tests on all commercial communications satellites. 
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4. Some test facilities have the capability to perform simultaneous multi-axis vibration testing, 
as needed. 

5. The vibration test provides qualification for tertiary/ancillary hardware that would not 
otherwise be tested. This includes: Cable harnesses, bellows, connectors, actuators, 
plumbing lines, wave guides, brackets, dampers, shades and shields, articulation & 
deployment mechanisms, shunt heaters, louvers, purge equipment, hinges and restraints, 
blankets & supports. 

6. The test provides an opportunity to determine the structural linearity in the operational 
vibration range of response. Note that linearity is typically assumed in coupled loads 
analysis (CLA). 

7. Force limiting reasonably accounts for the interaction with the base motion, and has been 
effectively employed in spacecraft vibration testing, thus reducing the potential for an over¬ 
test at the spacecraft’s natural frequencies in the test configuration. The force limiting 
takes into account the CI_A response levels. 

8. Force gauges under the spacecraft provide a very accurate method of measuring and 
limiting to the CLA loads during the vibration test for mid to high apparent mass modes. 

9. Numerous case histories have shown that vibration testing is effective for uncovering 
design or workmanship flaws which would have otherwise caused mission degradation or 
failure. 

21.4 Spacecraft Testing Experiences 

Note that NASA/GSFC typically uses sine vibration testing, whereas JPL tends to use random. 

A few examples from sine testing at GSFC are: 

• TRMM: During Observatory sine testing, found that the NASDA supplied PAF clamp band had 
insufficient tension and gapped during the test. As a result, the clamp band tension was 
increased for flight. 

• GOES had a workmanship problem involving a missing or loose bolt which caused structural 
failure of a mission-critical antenna. It was detected during the lateral sine test. 

• NOAA-K experienced IMU saturation during sine sweep testing. Because the spacecraft IMU 
provides guidance information for the Titan II launch vehicle during ascent, IMU saturation 
during launch would have resulted in a mission failure. Changes were made and launch vehicle 
restraints were implemented to resolve the problem, including wind restrictions at launch and a 
commanded first stage shutdown vs. fuel depletion. 

• TDRS-H: During the sine vibration test, the first two modes for the Space Ground Link antenna 
(SGL) were lower than predicted by the model. The first mode dropped from 15 Hz to 11 Hz 
and the second mode dropped from 33 Hz to 25 Hz. It turned out that the mathematical model 
of this “simple” antenna was wrong and therefore the Verification Loads Cycle had to be rerun. 

A few examples from random vibration testing at JPL are: 

• Cassini: Experienced an RTG electrical short to its spacecraft mount in system random vibration 
test. Significant degradation in spacecraft electrical power could have resulted. Spacecraft mount 
was redesigned. 
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• CloudSat: Cloud Profiling Radar waveguide failure in spacecraft random vibration test due to 
apparent poor workmanship of adhesive bonding. Possible loss of science data averted. 

• MER 1: Fundamental modes of the Rover in spacecraft random vibration test were 20% greater 
than predicted in all three axes. (Fixed base modal test had been performed on Rover, Lander, 
and Cruise Stage separately; FE models were then combined. Estimated stiffness of Lander 
attachment to Rover was too low.) FE model was updated just in time for the verification CLA 
cycle. Vibration test also revealed improper torque of bolts on some tanks in low level runs. 
Bolts were properly torqued and test completed successfully. 

• MSL Rover: experienced several motor encoder screws backed out of at least one of the Rover 
actuators during Rover random vibration test. The actuators are used throughout Rover and 
the issue was unlikely to have otherwise been found before launch, which could have been a 
serious threat to the mission. 

JPL prefers random vibration because it easier to control, particularly with respect to force limiting. 

Note that JPL tested the SMAP spacecraft to the following workmanship PSD: 

20 Hz to 250 Hz, 0.01 G A 2/Hz 

The overall level was 1.5 GRMS. The duration was one minute. The vibration was applied in the 
vertical axis and in one “45 degree” lateral axis. 

21.5 Sine Sweep Control 

Sine sweep vibration is more difficult to control than random especially for the case of lightly-damped 
modes. The sweep rate, compression factor, and tracking filter must be selected with great care. 

INVAP experienced control issues during sine vibration testing of the ARSAT-1 structural test 
model. But note that dynamic simulators were used for the test, which can have much less damping 
than the actual flight hardware. 
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22 Tacoma Narrows Bridge Failure 

22.1 Bridge Introduction 


Figure 22.1. Torsional Mode of the Tacoma Narrows Bridge 

The Tacoma Narrows Bridge failure is an example of self-excited vibration of such importance that it 
warrants its own section in this document. The original Tacoma Narrows Bridge was opened to traffic 
on July 1, 1940. It was located in Washington State, near Puget Sound. 

The Tacoma Narrows Bridge was the third-longest suspension bridge in the United States at the time, 
with a length of 5939 feet including approaches. Its two supporting towers were 425 feet high. The 
towers were 2800 feet apart. 

Prior to this time, most bridge designs were based on trusses, arches, and cantilevers to support heavy 
freight trains. Automobiles were obviously much lighter. Suspension bridges were both more elegant 
and economical than railway bridges. The suspension design thus became favored for automobile traffic. 
Unfortunately, engineers did not fully understand the forces acting upon bridges. Neither did they 
understand the response of the suspension bridge design to these poorly understood forces. 

Furthermore, the Tacoma Narrows Bridge was built with shallow plate girders instead of the deep 
stiffening trusses of railway bridges. Note that the wind can pass through trusses. Plate girders, on the 
other hand, present an obstacle to the wind. As a result of its design, the Tacoma Narrows Bridge 
experienced rolling undulations which were driven by the wind. It thus acquired the nickname 
"Galloping Gertie." 


An Introduction to Shock & Vibration Response Spectra: By Tom Irvine 


267 










enDAQ 


Chapter: Tacoma Narrows Bridge Failure 


22.2 Bridge Collapse 



Figure 22.2. Tacoma Narrows Bridge Collapse 

Strong winds caused the bridge to collapse on November 7, 1940. Initially, 35 mile per hour winds 
excited the bridge's transverse vibration mode, with an amplitude of 1.5 feet. This motion lasted 3 
hours. The wind then increased to 42 miles per hour. In addition, a support cable at mid-span snapped, 
resulting in an unbalanced loading condition. The bridge response thus changed to a 0.2 Hz torsional 
vibration mode, with an amplitude up to 28 feet. The torsional mode is shown in Figure 22.1. 

The torsional mode shape was such that the bridge was effectively divided into two halves. The two 
halves vibrated out-of-phase with one another. In other words, one half rotated clockwise, while the 
other rotated counter-clockwise. The two half spans then alternate polarities. 

One explanation of this is the "law of minimum energy." A suspension bridge may either twist as a 
whole or divide into half spans with opposite rotations. Nature prefers the two half-span option since 
this requires less wind energy. The dividing line between the two half-spans is called the "nodal line." 
Ideally, no rotation occurs along this line. The bridge collapsed during the excitation of this torsional 
mode. Specifically, a 600 foot length of the center span broke loose from the suspenders and fell a 
distance of 190 feet into the cold waters below. The failure is shown in Figure 22.2. 
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22.3 Failure Theories 

22.3.1 Candidates 

The fundamental weakness of the Tacoma Narrows Bridge was its extreme flexibility, both vertically and 
in torsion. This weakness was due to the shallowness of the stiffening girders and the narrowness of the 
roadway, relative to its span length. Engineers still debate the exact cause of its collapse, however. 
Three theories are: 

1. Random turbulence 

2. Periodic vortex shedding 

3. Aerodynamic instability (negative damping) 

These theories are taken from Reference [47]. Aerodynamic instability is the leading candidate. 

22.3.2 Random Turbulence 

An early theory was that the wind pressure simply excited the natural frequencies of the bridge. This 
condition is called "resonance." The problem with this theory is that resonance is a very precise 
phenomenon, requiring the driving force frequency to be at, or near, one of the system's natural 
frequencies in order to produce large oscillations. The turbulent wind pressure, however, would have 
varied randomly with time. Thus, turbulence would seem unlikely to have driven the observed steady 
oscillation of the bridge. 

22.3.3 Vortex Shedding 

Theodore von Karman, a famous aeronautical engineer, was convinced that vortex shedding drove the 
bridge oscillations. A diagram of vortex shedding around a spherical body is shown in Figure 22.3. Von 
Karman showed that blunt bodies such as bridge decks could also shed periodic vortices in their wakes. 

A problem with this theory is that the natural vortex shedding frequency was calculated to be 1 Hz, as 
shown in Section 22.4. This frequency is also called the "Strouhal frequency." The torsional mode 
frequency, however, was 0.2 Hz. This frequency was observed by Professor F. B. Farquharson, who 
witnessed the collapse of the bridge. The calculated vortex shedding frequency was five times higher 
than the torsional frequency. It was thus too high to have excited the torsional mode frequency. 

In addition to "von Karman" vortex shedding, a flutter-like pattern of vortices may have formed at a 
frequency coincident with the torsional oscillation mode. Whether these flutter vortices were a cause 
or an effect of the twisting motion is unclear. 
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Figure 22.3. Vortex Shedding around a Spherical Body 

22.3.4 Aerodynamic Instability 

Aerodynamic instability is a self-excited vibration. In this case, the alternating force that sustains the 
motion is created or controlled by the motion itself. The alternating force disappears when the motion 
disappears. This phenomenon is also modeled as free vibration with negative damping. 

Airfoil flutter and transmission line galloping are related examples of this instability. Further 
explanations of instability are given in References [48], [49] and [50]. 

The following scenario shows how aerodynamic instability may have caused the Tacoma Narrows Bridge 
to fail. For simplicity, consider the motion of only one span half. Assume that the wind direction was 
not perfectly horizontal, perhaps striking the bridge span from below, as shown in Figure 22.4. 



Figure 22.4. Span at Initial Rest Position 

Thus, the bridge is initially at an angle-of-attack with respect to the wind. Aerodynamic lift is generated 
because the pressure below the span is greater than the pressure above. This lift force effectively places 
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a torque, or moment, on the bridge. The span then begins to twist clockwise as show in Figure 22.5. 
Specifically, the windward edge rotates upward while the leeward edge rotates downward. 



Figure 22.5. Span Rotates Clockwise 



Figure 22.6. Span Rotates Counter-Clockwise 

The span has rotational stiffness, however. Thus, elastic strain energy builds up as the span rotates. 
Eventually, the stiffness moment overcomes the moment from the lift force. The span then reverses its 
course, now rotating counter-clockwise 

The span's angular momentum will not allow it to simply return to its initial rest position, however. The 
reason is that there is little or no energy dissipation mechanism. Thus, the span overshoots its initial 
rest position. In fact, it overshoots to the extent that the wind now strikes the span from above as 
shown in Figure 22.6. The wind's lift force now effectively places a counter-clockwise moment on the 
span. 
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Once again, strain energy builds up in the span material. Eventually, the stiffness moment exceeds the 
moment from the wind's lift force. The span thus reverse course, now rotating clockwise. Again, it 
overshoots its rest position. The cycle of oscillation begins anew from the position shown in Figure 
22.4, except that the span now has rotational velocity as it passes through the original rest position. 

The cycles of oscillation continue in a repetitive manner. 

Note that the wind force varies as a function of the span angle during the cycle. The wind force may 
also vary with the angular velocity. The wind force is not a function of time, however. 

Eventually, one or two failure modes occurs. One possibility is that the span experiences fatigue failure 
due to an excessive number of stress reversals. The other is that the angular displacement increased in 
an unstable manner until the material is stressed beyond its yield point, and then beyond its ultimate 
stress limit. 

These two failure modes can be interrelated. For example, accumulated fatigue effectively lowers the 
yield and ultimate stress limits. Regardless, the bridge collapses. 

As a final note, the aerodynamic instability oscillation is not a resonant oscillation since the wind does 
not have a forcing frequency at, or near, the bridge's torsional mode frequency. Some physics and 
engineering textbooks mistakenly cite the Tacoma Narrows Bridge as an example of resonance. This 
problem is discussed in Reference [51]. Nevertheless, the bridge's collapse remains the most famous 
structural failure due to vibration. 

22.4 Bridge's Strouhal Frequency 

The original Tacoma Narrows Bridge collapsed in 1940. It experienced severe torsional oscillations 
driven by a 42 mile per hour wind. The fundamental weakness of the Tacoma Narrows Bridge was its 
extreme flexibility, both vertically and in torsion. This weakness was due to the shallowness of the 
stiffening girders and the narrowness of the roadway, relative to its span length. A cross-section of the 
bridge is shown in the next figure. 


39 ft 



Figure 22.7. Cross-section of Tacoma Narrows Bridge Span 

The Strouhal number for the bridge cross-section is S = 0.11 according to Reference [51]. 
Furthermore, the characteristic dimension is the girder height 8 ft, which may be substituted for the 
diameter. The Strouhal frequency for a 42 mph (61.6 ft/sec) wind is thus 

f s =S U/D (22.1) 
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f s =(0.11) (61.1 ft /sec) / (8 ft) =0.84 Hz (22.2) 

The Strouhal frequency is rounded to 1 Hz in Reference [51]. On the other hand, the bridge's torsional 
oscillation was observed to be 0.2 Hz. The Strouhal frequency was thus at least two octaves greater 
than the torsional mode frequency. Billah and Scanlan thus argue in Reference [51] that the Karman 
vortex shedding could not have driven the torsional mode oscillation. They argue instead that an 
aerodynamic instability resulted in a self-excited oscillation, which caused the failure. Several decades 
after the failure, the matter is not completely settled. 

22.5 Replacement Bridge 



Figure 22.8. New Tacoma Narrows Bridge 

A new Tacoma Narrows Bridge was built in 1950, as shown in Figure 22.8. The new bridge had truss- 
girders which allowed the winds to pass through. It also had increased torsional stiffness because it was 
thicker and wider. Furthermore, wind tunnel testing was performed to verify the design of the new 
bridge prior to its construction. A similar bridge was later built alongside the new bridge in the figure 
to ease traffic flow, with one-way traffic on each bridge. 
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